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Engineers, anxious to base their calculations on a rigorous as well as simple theory, 
do not hesitate to assume that structural materials are perfectly elastic, that is, that 
the deformations experienced under the action of external loads disappear as soon as 
these loads are removed. But this assumption is justified only if the stresses remain 
everywhere below a certain limit, called the “elastic limit’’, and if the period over which 
the loads are applied is sufficiently short. 

In practice, the stresses exceed the elastic limit much more frequently than is gener- 
ally supposed, even in the most rigorously and safely calculated structures. The duration 
of the stresses, on the other hand, is extremely variable, covering all the range from 
impact stresses to stresses which subsist during the entire lifetime of the structure. 

Stresses beyond the elastic limit, as well as long sustained stresses, produce non- 
elastic strains which do not disappear on removal of the loads which have caused them. 
This fact has two important consequences. First, the assumption of a one-to-one corre- 


spondence between strains and stresses, which is fundamental in the theory of elasticity, 
must be abandoned. Secondly, a state of residual stresses** is established, causing a 
distribution of the inner stresses which is completely different from the distribution 


furnished by the theory of elasticity. 

Obviously, we cannot represent this phenomenon by a mathematical theory without 
introducing some hypothesis regarding the nature of the deformation and its relations 
to the stresses which produce it. This is what we are proposing to do next, and we will 
try to find the physical significance of these hypotheses and point out to what extent 
they correspond to the empirical facts. 

A solid body is said to be perfectly elastic if there exists in it a one-to-one corre- 
spondence between stresses and strains. This is the case for most structural materials, 
as long as the stresses are sufficiently small and act for a sufficiently short time. But, 
as soon as the stresses become greater, though they may act over a very short period, 
the process ceases to be reversible; part of the strain will remain in the structure even 
after the loads have been removed. We call these “plastic strains’’. 

Of course, we do not intend to assert that there actually exists, for each material, 
a well-determined limit such that the material behaves in a perfectly elastic manner 


*Received April 4, 1949. This paper was presented at the Second Symposium on Plasticity sponsored 
jointly by the Office of Naval Research and the Bureau of Ships (Brown University, April 4-6, 1949). 

**The term “residual stresses” as used in the following refers to the state of stress as well as to 
the state of permanent strain which causes the residual stresses. It would be desirable to adopt a single 
English term for this concept corresponding to the Italian coazione or the French coaction. 
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as leng as the stresses do not exceed this limit. We wish to state only that, for all practical 
purposes, we may assign a limit, below which the process can be treated as reversible, 
in the sense that the plastic strains are negligible in comparison with the elastic strains. 
Beyond this limit, however, the process becomes decidedly irreversible. In general, this 
limit is purely conventional. As a matter of fact, the limit, beyond which irreversible 
processes are observed, decreases with increasing accuracy of the observations. 

However, there are materials, such as iron, mild steel, bronze, etc., for which the 
limit between the elastic range (characterized by a nearly complete absence of permanent 
deformation) and plastic range (where almost the whole deformation is permanent) is 
well-defined. In these materials the permanent strains, though occurring also when 
the stresses are very small, are negligible compared with the elastic strains as long as 
the stresses do not exceed a certain critical value. As soon as this value is reached, the 
laws which govern the process change all of a sudden, and we may say that all further 
strains are practically permanent. As a matter of fact, the plastic strains which now 
occur, grow so rapidly, even under a constant stress, that the elastic strains which 
affected the material up to this moment, become negligible in comparison. This limit, 
which is no longer merely conventional but has a well-determined physical meaning, is 
called the “yield limit’. 

We will choose this apparently rather special case as the starting point to derive 
a mathematical theory of the elasto-plastic equilibrium, that is of the equilibrium of 
bodies in which the elastic limit has been reached and plastic deformations have oc- 
curred. 

We will make the following assumptions: first, the one-to-one correspondence be- 
tween stresses and strains is valid for each element of the body as long as the material 
is stressed within its elastic range; secondly, as soon as this limit has been reached, 
the element becomes perfectly plastic and the deformation can increase indefinitely 
under constant stresses. Moreover, we will assume that the plastic strains are super- 
imposed on the elastic strains without changing the characteristics of the latter, partic- 
ularly without depriving them of the property of disappearing upon the removal of the 
loads which have caused them. 

It can be shown that these apparently very restrictive hypotheses adequately cover 
a great number of cases. They even explain the cases in which the transition from the 
elastic zone to the plastic zone takes place with continuity, that is, when a physically 
well-determined elastic limit is absent. From these hypotheses it is also possible to 
develop a very elegant interpretation of the so-called “strain-hardening’’, which is 
always connected with plastic deformation. 

To this end we will first investigate the nature and the fundamental characteristics 
of this phenomenon, referring to a very particular, but simple and experimentally well- 
known case: that of a cylindrical or prismatical metallic bar submitted to axial tension. 
It is well-known that when the yield limit is reached in such a bar, localized slip occurs 
along certain planes, whose angle with the axis of the bar is determined by a relation 
between the normal and the shearing components of the stress transmitted across these 
planes. 

A single, small inhomogeneity in the material is sufficient to start localized slip. 
In the zones where slip occurs it has the tendency to cause a reduction of the cross 
section together with an increased strength of the material. At first, the favourable 
effect of this second factor more than balances the unfavourable consequences of the 
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first one. Therefore, slip soon ceases in that section only to appear immediately in some 
other section, then in a third one, and so on. This process of propagation may be made 
visible by polishing carefully the surface of the test specimen. Slip is then indicated by 
the appearance on the polished surface of extremely thin inclined lines, the so-called 
slip bands. These bands grow in number, extend, and become interlaced in a kind of 
opaque veil which quickly spreads and soon covers the whole surface uniformly. 

We are thus led to consider a new condition of the material in which the original 
homogeneity, obviously destroyed by localized slip, is finally restored on the whole. 
The most remarkable and obvious characteristic of this new condition is a rise of the 
elastic limit of the material, that is, an increase of the value of the stresses to which 
the test specimen can be submitted without showing plastic strains. It is exactly this 
rise of the elastic limit which is called the strain-hardening of the material. 

Having reached these conclusions from the experimental standpoint, we are now 
going to look at the problem from an exclusively theoretical point of view. We will 
assume that the material is perfectly homogeneous and consequently that in the elastic 
field the stresses are uniformly distributed throughout the specimen. According to this 
hypothesis, the elastic limit will be reached in the same instant in all parts of the body, 
and we may conclude that the plastic strains will be uniformly distributed too. 

These plastic strains will not give rise to any residual stresses; they will be super- 
imposed on the elastic strains, without changing the stress distribution or the value of 
the strain energy. This means that the mechanical work necessary to produce these 
plastic strains has been completely spent to overcome the inner friction and to produce 
heat. 

Let us now suppose that the specimen is unloaded. The conclusions reached above 
explain the fact that the specimen has a permanent set, in the strictly geometrical 
sense of a permanent change of shape and dimensions, but they do not indicate at all 
the experimentally established raising of the elastic limit. 

However, due consideration must be given to the fact that the assumed structural 
homogeneity of the material, as well as the uniform distribution of the plastic strains, 
do not really occur. It then becomes clear that, since the plastic strains are no longer 
“compatible” as a rule, a state of residual stresses is set up, which will superimpose 
itself on the stresses previously established in the elastic range and which will subsist 
even after all elastic strains have disappeared. In order that this condition may be 
achieved, at least a portion of the mechanical work spent on producing the elastic strains 
should not be lost, but transformed into potential energy. This energy will not be re- 
covered as external work when the loads are removed. For this reason we shall call it 
“residual energy”’. 

However, the simple fact that this residual energy subsists after the external forces 
have ceased acting, is sufficient to let us think that something else beyond shape and 
dimension has changed in the body, and that here lies, perhaps, the explanation of the 
phenomenon with which we are dealing. We will show that this is possible by means of 
a very simple but very significant example, taken from the field of reinforced concrete 
design. 

Let us consider a reinforced concrete bar submitted to an axial tension. It is well- 
known that, as long as the stresses are small and the material behaves everywhere in a 
perfectly elastic manner, the stresses must be distributed over the various elements of 
each cross section of the bar according to the ratios of their moduli of elasticity, if the 
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cross section is to remain plane. But as soon as the elastic limit of concrete is reached, 
the material enters the plastic state. Again, if the cross section is to remain plane, the 
concrete must undergo, independently from any further increase of the stresses, such 
strains as are necessary to let the metallic reinforcement carry any further increase of 
the load. 

The stress-strain diagram then takes the shape of a broken line, of which the first 
section, beginning at the origin, covers the elastic range. The second section covers the 
elastic-plastic range. A third section, parallel to the strains axis will then represent the 
final phase, when thc entire cross section behaves in a plastic manner. 

Let us now suppose that the load is removed. The elastic strains will naturally show 
a tendency to disappear; this will not be the case for plastic strains. On the contrary, 
these will even succeed in preventing a complete disappearance of the elastic strains. 
Thus, a state of residual stress will be set up, such that the metallic reinforcement will 
be stressed in tension, while the concrete will be under compressive stresses. It follows 
from this that the load required in order to let the concrete again reach its elastic limit, 
so that further plastic strains occur, must now be greater than the first time. 

If we ignore the existence of the residual stresses, we shall naturally speak of a rise 
of the tensile elastic limit of the concrete. But this is not correct: it is not the elastic 

, 


limit of the concrete which has risen, it is the initial state of stress which is no longe 
the same. 

This result is not at all restricted to the type of structure we have considered as 
an example. To make sure of that, let us go back to the fundamental problem, from 
which we started, and study the general case of any structure whatsoever, consisting 
of heterogeneous elements having elastic limits varying from a certain minimum to 
some maximum. Assuming these elements to be uniformly distributed throughout the 
structure, we can see in a very general way that the above-mentioned broken line will 
be replaced by a polygonal line and, as a limit, by a curve. In fact, each element stops 
contributing to the resistance of the structure as soon as it reaches its own elastic limit, 
and thereafter it just conforms plastically to the elastic strains of the surrounding 
elements, which have higher elastic limits. These elements must carry any further 
increase of the load. In other terms, the stress is shifted from the elements which have 
reached their yield limit to the elements which are still stressed elastically. If now the 
loads are removed, the plastic strains prevent the structure from returning to the 
original unstrained shape and thus give rise to residual stresses. The elements which, 
owing to their higher yield limit, have carried the greatest portion of the load, keep 
on being subject to stresses of the same type, whereas the elements which have under- 
gone plastic strains show a tendency to be subject to stresses of the opposite type. 

It follows from this that when the structure is again loaded, the last elements are in 
a more favourable condition, because the residual stresses will have to be overcome 
before stresses of the type connected with the acting loads may again appear. Thus the 
elastic limit shall be reached, even for the weaker elements, when the loads have attained 
a higher value than was originally necessary to produce plastic strain. This explains the 
apparent rise of the elastic limit of the material. 

If, on the other hand, after having produced the strain hardening of a specimen, 
by submitting it to the action of a well-determined system of loads, we apply to the same 
specimen a system of loads of the opposite direction, the stresses produced by these 
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forces will be added to, and not subtracted from, the residual stresses. It will then 
appear that the elastic limit has been lowered (Bauschinger effect). 

We are thus entitled to consider our interpretation of the phenomenon of strain- 
hardening as reasonable. At any rate, it is certain that the outline of the phenomenon, 
which we have suggested as a starting point of a mathematical theory of the elasto- 
plastic equilibrium, leads to a very satisfactory interpretation of the experimental 
results. 

The situation is very different when the stresses, though not being high enough to 
cause the elastic limit to be reached in any point of the structure, keep on acting for a 
time long enough to result in irreversible effects, that is to give rise to strains which 
are not going to disappear immediately when the loads producing them are removed. 
This phenomenon has been long known as ‘‘time effect’? and involves the superposition 
on the elastic strains, due to the external forces, of other strains which grow progressively 
with time under the sustained action of those same forces. These strains do not disappear 
as soon as the forces cease acting, but they may disappear later on (in any case very 
slowly) if the structure remains free to return to its initial unstrained shape. 

A mathematical theory of such a phenomenon is not easy to establish, because the 
strains at a certain instant depend not only on the loads acting at this instant, but also 
on previous loads and on the duration of the period for which they have acted. We 
shall consider only a very particular case, which is, however, very important for practical 
purposes, namely, the case of very small strains (of the magnitude of the elastic strains) 
which we shall assume to increase, at any instant, in proportion with the stress existing 
at this instant. 

Let us suppose that the structure is homogeneous and that the coefficient of pro- 
portionality, which is usually called the “creep coefficient’’, is the same for the entire 
structure. Moreover, let us assume that there is no state of initial stress. The time effect 
will then involve strains which, like the elastic strains, satisfy the conditions of com- 
patibility. Accordingly, these strains will subsist even when the elastic strains will have 
disappeared with the removal of the loads. It is clear, therefore, that the time effect 
does not produce any residual stresses, if the structure was free of initial stresses. It 
can be shown quite generally that in a homogeneous body the time effect appears geo- 
metrically as a proportional increase of the elastic deformation produced by the loads 
and statically as a proportional decrease of the residual stresses. 

Let us now introduce the hypothesis that the body is not homogeneous, which 
proved so useful when we dealt with plastic strains. The situation will then change, and 
the theory will explain immediately the observed tendency of the strains, due to time 
effect, to vanish slowly once the loads have been removed. 

We will assume, therefore, that the body consists of heterogeneous elements affected 
by different flow coefficients, uniformly distributed throughout the body. We may 
foresee that the elements characterized by higher flow coefficients, following a prolonged 
action of the external forces, will undergo greater strains than the elements with lower 
values of that coefficient. As a consequence, the last elements, conforming to the first, 
will undergo supplementary elastic strains. When the loads are removed, these addi- 
tional strains, like every elastic deformation, will tend to vanish, but they will be pre- 
vented from vanishing completely by the non-elastie strains of the surrounding elements. 
This gives rise to a state of residual stress: the stresses caused by the elements with 
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small non-elastic strains will produce in the elements with greater inelastic strains a 
new set of strains of opposite type, and these will force the body to move towards the 
initial unstrained shape. 

From this mass of apparently so different and complicated facts a very interesting 
and absolutely general conclusion can be drawn concerning the co-existence, at any 
moment and in any state of equilibrium, of elastic strains and non-elastic strains, what- 
ever their origin and their physical nature may be. Here is how we can reach this con- 
clusion quickly. 

We will denote the components of the elastic strain by 

eg Sy "°° 9 Tp 


and those of the plastic strain by 


We will further assume that both strain distributions are continuous, uniform and very 
small, so that it is always possible to write the components of the resulting total strain 
as follows: 
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It must be kept in mind that the resulting total strains are compatible, but that 
neither one of the two partial strains will be compatible as a rule. Assuming that neither 
one of the two partial strains is compatible, let us now try to state analytically that one 
of these partial strains is elastic. 

This can be achieved in two ways, either considering the deformation as a state of 
equilibrium, or taking into account the way it occurs or vanishes. 

In the first case we are led to assume that any element of the body is kept in a 
strained condition by the stresses which act on it owing to the surrounding elements. 
Thus the strains are strictly connected with a system of stresses which we shall define, 
as usual, by means of six components of stress, namely: 

gg OFF oh Tees 


If, on the other hand, we adopt the second point of view, taking into account the 
tendency of the strains to vanish once the loads which caused them have been removed, 
we must conclude that there exists an elastic strain energy, produced by the trans- 
formation of the mechanical work spent to make the actual strains occur in the body. 
This energy may again be recovered as external work, when the strains disappear for 
one reason or the other. This strain energy is assumed to be the sum of the elementary 
strain energies of all the individual elements of the body. In this way each of these 
elementary energies is completely defined when the condition of the element to which 
it refers is known. It is not necessary, in order to characterize it, to introduce either the 
condition of the other elements or their relative position. This results in the following 
expression of the energy of the body: 


a 
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where ¢ is the elastic strain energy referred to the unit of volume. By definition, ¢ is a 
function of the six strain components. 

It can then be shown that the elementary elastic strain energy is a positive definite 

form of the second order, whose first partial derivatives with respect to the six strain 

components equal the corresponding stress components. Thus the assumption from 
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which we started amounts to this: the six stress components are linear and homogeneous 
functions of the six strain components, and vice versa. It follows from this that the 
elementary elastic strain energy may also be regarded as a second-order, positive definite 
form of the six stress components. Its first partial derivatives with respect to these new 
variables equal the corresponding strain components. 

Let us now study the non-elastic strains, that is, the strains which do not vanish 
together with the loads. The work spent in order to produce these strains will not be 
recovered immediately. To develop a better understanding of this we will study sepa- 
rately the element, taken alone, and the body consisting of all elements. 

In the first case, the fact that the mechanical work spent to produce the strains 
cannot be recovered may be expressed by assuming that this work is lost through ir- 
reversible transformations. This prevents the accumulation of any strain energy. The 
same idea can equally be expressed by stating that the strains, once they have been 
produced, go on existing independently from any external action which may keep the 
body in its strained condition. As a consequence, the existence of strains does not imply 
the presence of stresses. 

But what has been stated for a single element does not apply to the body as a whole. 
In it non-elastic strains can occur only in some points, so that they will generally form 
a non-compatible system. Along with this system we must then have a series of com- 
plementary elastic strains such that the composite strains become compatible. Thus 
the non-elastic strain of a given element gives rise to a series of elastic strains and of 
stresses in the surrounding elements. We shall call these stresses “‘residual stresses’’ to 
distinguish them from the stresses which are connected with the loads. These stresses 
involve an elastic energy which we shall call “residual energy”’ to remind us that it is 
connected, like the condition to which it refers, with the presence of non-elastic strains 
and that it keeps on existing even after the loads which produced it. 

The mechanical work required by a non-elastic deformation does then consist of 
two parts: one is spent to produce the deformation itself and is lost through irreversible 
processes; the other gives rise to the state of residual stress. Neither of these two func- 
tions is recovered the moment the loads are removed. But while the first one is definitely 
lost, the second will remain available in the body as strain energy. It may even happen 
that, joining the strain energy produced in the same body by another system of loads, 
this residual energy will become apparent through the reactions of the body. To the 
new loads, it may even be transformed again into mechanical work, if there exist in 
the body any elements able to be further strained under time effect, or if by proper 
treatment, we succeed in freeing it. In any case, this distinction between the two parts 
of the work spent on the production of non-elastic strains is of basic importance. It 
will always be necessary to take into account separately the part which causes a real 
loss of energy, while, in view of global evaluation of the strain energy, the other fraction 
will have to be treated as if it were work spent in order to produce new elastic strains. 

Let us now consider the body in equilibrium under the action of a given system of 
loads and suppose that the stresses undergo a very small variation compatible with 
these loads. This means that, while the non-elastic strains remain unaltered, the stress 
components shall undergo very small variations 
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which form a system of stresses in equilibrium in the absence of external forces. 
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We consider the function 
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the first variation of which can be written in the form: 
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We see at once that this variation yields the work which the system of stresses 
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would carry out if the body were to undergo the change of shape characterized by the 
following components: 
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Now this system of stresses is, by hypothesis, in equilibrium in the absence of all external 
forces. On the other hand, the change of shape is certainly compatible with the possible 
constraints of the body, since it is exactly the one which the body had to undergo in 
order to pass from the unstrained natural state on to the state of equilibrium which we 
are considering. From the principle of virtual works it follows, therefore, that 


Taking into account the fact that the second variation is: 
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and that the integrand is positive definite, we may state the following theorem: For 
each system of loads and of non-elastic strains, the inner stresses characterizing the state of 
equilibrium are those which minimize the expression 
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If there are no non-elastic strains, the problem fits into the framework of the classical 
theory of elasticity, and the function which has to be made a minimum is the elastic 
strain energy. We find here again, as it could have been foreseen, the theorem of Mena- 
brea. 

It is of basic importance to note that we can pass from the equation which is the 
basis of the classical theory of elasticity, to the equation which appears to be a possible 
basis of a new theory of the elasto-plastic equilibrium, simply by replacing the six 
elastic strain components by the six components of the total strain. Provided an analogy 
exists between the surface conditions involved in each particular problem, the solution 
to which the classical theory of elasticity leads for its problems may, therefore, be 
transformed into as many solutions of analogous problems of elasto-plastic equilibrium. 

In this way we are able to establish a new theory of the strength of beams, the most 
important conclusion of which is a new outlook concerning the redundancy of constraints 
and the consequent statical indeterminacy of the structures. By the classical theory of 
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elasticity one has become accustomed to consider statically indeterminate structures as 
extremely sensible to all non-elastic strains. This sensibility was often considered as a 
disadvantage, sometimes even as a real danger, since non-elastic strains were beyond 
the scope of our calculations and, consequently, of our forecasts. Today, as this scope 
has extended to cover the whole field of non-elastic strains and we are able to take exact 
account of the effects of these strains and to calculate the changes they will cause in 
the state of equilibrium of the structures, we should consider this sensibility rather as 
an advantageous and even precious quality of the hyperstatic structures. 

Let us consider, for instance, the case of proper plastic strains, that is, the strains 
which occur in the outermost fibers of a beam as soon as the elastic limit of the material 
has been reached. Our theory points out in the most simple and clear way three types 
of effects which such strains may have. The first type concerns exclusively the distribu- 
tion of stresses over the section itself and becomes apparent by limiting the stresses in 
the fibers where the elastic limit has been reached: this limitation is compensated for 
by an increase of stresses in the inner zone of the section, where the strength of the 
material was less utilized originally. This new stress distribution is due to the static 
indeterminacy of the cross-section itself; it takes place, therefore, even when the beam 
is statically determinate, as far as external constraints are concerned. The second type, 
on the other hand, involves the static indeterminacy of the diagram of bending moments, 
that is, the existence of at least one redundant constraint; it becomes apparent by 
limiting the value of the moment in the section where the plastic strains have occurred 
and by increasing it in the sections which were less stressed originally. The third type, 
finally, is observed only if the constraints involve further static indeterminacy of the 
structure; it becomes apparent through a force which modifies the distribution of the 
stresses not in one or more se¢tions, where plastic strains have occurred, but in all the 
sections of the beam. This force tends to limit that, of the two maximum stresses, which 
first has reached the elastic limit, and to increase the other. 

It is by a proper combination of these effects, that is, by using each of them in the 
most suitable way and by following laws which only a rigorous and complete analysis 
of the phenomenon can discover, that nature carries out, within the limits fixed by the 
data of the problem, the best possible utilization of the resistance of the material. Be- 
sides, it is clear that this behaviour of the structure, which is now at our disposal since 
we know its secret, will have a wider field of application and will allow us to obtain 
important simplifications of the analysis of stresses in indeterminate structures. There- 
fore, we must get used to looking at the problem in a new way. 

The classical theory had led us to think of the state of equilibrium of a system sub- 
mitted to given loads as something well-determined in terms of these loads. Even in 
the statically indeterminate case, in which the equation of statics yields an infinite 
number of solutions, the hypothesis of perfect elasticity of the material is sufficient to 
determine the solution completely. This result does not represent the truth, however. 
Whenever elastic strains are accompanied by even very small non-elastic strains, it is 
no longer true that only one solution of the equations of statics is physically possible. 
On the contrary, they all are possible, and which one shall come into being depends 
upon the occurring of suitable non-elastic strains and of the state of residual stress 
connected with them. 

But here the field of application widens. The question is no longer about the plastic 
strains only and about the part they can play in determining the state of equilibrium. 
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As a matter of fact, the study of plastic strains has allowed us to catch a glimpse of 
possibilities which can be accomplished in quite a different way. It is sufficient to recall 
the modern technique of the prestressed materials, especially of prestressed concrete, to 
get an idea of the wide applications of the theorem which we have just derived and of 
its use in the study of the strength of materials. 


Acknowledgement: The author thanks Ing. Pier Luigi Romita, who has prepared the 
English text of this paper, 
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THE B.W.K. APPROXIMATION AND HILL’S EQUATION, II.* 
BY 


L. BRILLOUIN** 
Cruft Laboratory, Harvard University 


1. The B.W.K. Method. The B.W.K. procedure’ was discovered in connection with 
problems of wave mechanics, where Planck’s constant h was treated as a very small 
quantity, and functions could be expanded in series of powers of h. The first terms of 
the expansion correspond to classical mechanics (completed with some quantum condi- 
tions) and higher terms represent typical wave-mechanical effects. Optical problems can 
be discussed along similar lines, starting with “geometrical optics,” and obtaining 
‘“nhysical optics” (diffraction, for instance) as higher order corrections. 

The B.W.K. method was discussed recently” by the present author in connection 
with its possible use in a purely mathematical problem, the discussion of Hill’s equation. 
The method sketched on this occasion proved very valuable and showed the need for 
a more complete discussion of the whole question. It will be shown that the B.W.K. 
procedure can yield a very good approximation in a great many mathematical problems, 
and leads directly to asymptotic expressions similar to those obtained by P. Debye in 
the case of Bessel functions. 

There is a variety of problems leading to equations of the general type 


da* 
where F(x) is a real function. We may, for instance, discuss wave propagation along a 
linet whose properties vary from place to place. If W(x) is the velocity of propagation 
and w the frequency, the solution is 
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+The propagation of waves along a dissipationless electric line leads to the following equations 
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where in general /(r) and c(x) are two functions of x. One may introduce a new variable y defined by 
l, dy = l dz, (A.2) 
where Jp is a constant that may represent the average value of I(x). This yields 
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an equation which corresponds to our general type (1). It should be noticed that (//c)!/ represents the 
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fe'*' (2) 


with 
“ar 


P= 


f satisfying Eq. (1). This representation will be found very useful and will enable us 
to give a simple physical translation of our mathematical formulas. 
Let us now consider a function 


ax 


w= GV, S= | Gaz (3) 


“ 
a 


where G(x) is any given function of x. Straightforward computation shows that u 


satisfies the following equation: 


ws , sfeot.rte° 
u +uf 6 - 5 (=) +12] = 0, (4) 


and a comparison between (4) and (1) leads to the condition 
» Bey , 1e" 
F=G4 - + i (5) 
4\G 2G 
If this equation can be solved to a certain approximation, the function u of Eq. (3) 
represents an approximate solution of Eq. (1). On many occasions, terms in G’ and G’”’ 
can be omitted in Eq. (5), and a reasonably good approximation is obtained by taking 


G=G, = F™. (6) 


This is the zero order approximation in the B.W.K. procedure, and the question is to 
discuss its limits of validity. 

A special case of importance was presented in the preceding paper (L.B.H. 1, p. 172). 
With 
_A -_ 
(a+ 2x)’ 


‘ cancel out in Eq. (5) making (3) with (6) a rigorous solution of Eq. (1). 


G = (7) 


terms in G’, G’ 


characteristic impedance of the line, hence the result: if /(x) and c(x) vary by the same amount, and c/ / 
remains a constant, Eq. (A.4) reduces to the standard type of wave equation. There is no reflection and 


no perturbation to wave propagation. 
The change of variables (A.2) means that we measure the length of the line in “henrys”’ and not in 


yards. We could also measure the length in “farads,” by taking 

c, dz =cdz (A.5) 
and thus obtain another equation: 
a] 2 Lal = 
Dz" ae 


When / and c vary in such a way as to keep a constant ratio (constant characteristic impedance) Eq. (A.6/ 


0 (A.6) 


Oz 


is again a standard wave equation, and the variables y and z are proportional to each other. 
These remarks show the connection between our discussion and the one given by Schelkunoff (5) 
£ : 


who considers the case of two independent functions /(x) and c(z). 
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2. Systematic successive approximations. A systematic set of successive approxima- 
tions can be worked out under the following assumptions: 


|G | ay |). , 

| 72 | < & ~s i <€, e <1. (8) 
| Go | G5 | 
We may look for an expansion 

G= Go + G, oe G, + --- (9) 


where the successive terms would be of order e’, e’, --- , and substitute into Eq. (5): 


1/2 
G, =F”, 


26; 3 (Ge) L= 





“4\G,/ ~ 2G,’ 
00.6, + @ = 3(G&)(G _ G) _ 1G (Gi _ Gs, 
a ee (#) (@ ~ “) -— 3G (& -G) (10) 


The expansion G of Eq. (9) would then yield, according to Eq. (3), a first solution of 
Eq. (1), 
u= Ge S = [ Gdz. (3) 


The next step is to obtain another independent solution v of Eq. (1). Two cases 
must be distinguished here: in the first case, 
(A) F > 0, G real and positive, S real; (11) 


y-1/2_+:i8 
v=yu*=G'e’’, 


the complex conjugate of u, gives a second independent solution. If we think of a problem 
of propagating waves (Kq. 2), wu corresponds to a wave propagating to the right, and v 
to a wave propagating to the left. The second case arises when 
(B) F < 0, Gand S purely imaginary, 7'’*u real. (12) 
In that case, we may use the general relation (12) of the preceding paper (L.B.H. 1, 
Eq. 12), 
uv’ -- vu’ = C, (13) 

which must hold between two independent solutions of Eq. (1). The constant C can be 
taken as unity for convenient normalization.* When wu is given, Eq. (13) represents a 
first order differential equation for v and its solution is easily found as 

e 

ad . 
v= Cu [ “ de. (14) 
“0 

The first case (11) corresponds to propagating waves; the second one (12) represents 


attenuated waves without propagation. 


*In the preceding example, where v = u* the constant C equals 27. 
I £ 
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In any region x where our approximate solutions u, v can be obtained, the general 
solution of Eq. (1) reads 


f = Au + Bu. (15) 


The difficulty begins with the problem of connecting these different regions and choosing 
in each of them coefficients A, B which will correspond to approximations to one and 
the same solution of Eq. (1). This problem was first clearly stated by Kramers and 
discussed by Kemble, Furry and other authors. We are going to discuss that problem 
again and show that many important cases have been overlooked by previous authors. 
3. Cases of exception—Connection of different regions. The method developed in 
Sec. 2 is based upon assumptions (8) and stops working in the following three cases: 


. 


G, — 0, (16) 
Gi! +00. (18) 


We must now discuss what happens when such critical conditions are realized. 
Case (16) is apparently the only one which has been discussed by previous authors. 
In the neighborhood of a point where the function F (Eq. 1) is zero, one may expand 


F = b(4 — %)", n>0 (19) 


and discuss the corresponding equation. Taking x, = 0 brings the zero point to the origin, 
and Eq. (1) reduces to 


f’ + bf =0 (20) 


whose solution is (Jahnke-Emde tables, Dover ed. p. 147) 


j\/2r 2(b)'”* (n+2)/2 ‘ 
f = 2 Zi, n+2) (? + 2 ’ (21) 





where Z is a solution of Bessel’s equation of order 1/(n + 2). This covers the case 
originally discussed by Kramers, when the F function has a single root at the origin 





f= 2,20 2°”), (22) 
Assuming b > 0, we consider both sides of the origin: 

(a) z>0, F>0, 4 = Hizor Hin. (23) 
Here the functions corresponding to the two B.W.K. solutions u, and u¥ (Eq. 11) are 
the two Hankel functions. On the other side of the origin, 

(b) -<@ F <0, Zia = Jig OF Yiy3 . (24) 


Here the functions corresponding to the B.W.K. solution u,; and to the second solution 
v, (Eq. 14) are the Bessel and Neumann functions. 
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Hence the B.W.K. approximate solution u, to the right is not the analytic con- 
tinuation of the B.W.K. solution u, to the left, but 


u, goes over into 1/2 (u, + u*). (25) 


A more complete discussion of this case is found in Kemble’s book’ and in Furry’s 
paper. A refinement of that discussion was recently given by Isao Imai in a very in- 
teresting letter to the Phys. Rev. 74, 113 (1948). If we take the case when Eq. (1) results 
from a problem of wave propagation, as stated in Eq. (2), condition (16) corresponds 
to the total reflection at the end of the line. This is apparent since u, and u* represent 
two waves propagating in opposite directions, the superposition of which is needed to 
match the u, attenuated wave in the region where propagation is no longer possible. 
This point of total reflection is often called a “turning point.”’ 

Cases (17) and (18) correspond to rather sudden changes in the properties of the 
line, with propagation on both sides of the discontinuity (provided F(x) remains positive). 
This is a case of mismatch at a junction and results in partial reflection as we shall see. 


F 








0 
FIG. 1 


Let us consider first the case of condition (17), assuming that F’ becomes very large 
at x = 0. This corresponds to Fig. 1 with an inflection point at K. Simplified examples 
are shown in Figs. 2 and 3, when the function F has a discontinuity 2H at the origin, 

SLOPE b 


K+H K+H 
dt il 


| 2H 




















1 
{ 
cul mi 
ee proadomSe 
“ i 
ie) Xo-a O +a « 


16. & FIG. 3 


or a very rapid change from K — H to K + H on a small interval —a < x < a. In 
the case of Fig. 2, with a discontinuity in F, we have a solution on the left 


fi = Au, + By, (26.a) 








368 L. BRILLOUIN [Vol. VII, No. 4 


to match a similar solution on the right 
f, = A,u, + B,v, . (26.b) 
We must write the continuity of f and f’, and this determines A, , B, when A, , B, are 


given, or vice versa. This was the procedure followed in the preceding paper (L.B.H. 1, 
p.169). Assuming the function F(x) to remain positive on both sides of the origin, 


u and v = u* correspond to waves traveling in opposite directions. If we have a single 
wave on the right (B, = 0) we need two waves on the left. Hence an incident wave 


from the left generates a reflected wave traveling backwards and a transmitted wave 
penetrating into the region to the right. This obviously means partial reflection at the 
junction z = 0. 

In the example of Fig. 3, we would use a solution (26.a) for x < —a and (26.b) 
for x > a, together with solutions of the type (22) in the interval —a < x < a: 





1/2 (1) (2 : 2(b)'”” 3/2 - 
f =2”[AH;,3(y) + BHy3(y)], aioe (x — 2%)”. (27) 
We must write two continuity conditions at « = —a and two more at x = a, thus 


enabling us to compute the 4 coefficients A, B, A, , B, when A, , B, are given. The 
net result is again partial reflection, and this is what should be expected to happen in 
the more difficult problem of Fig. 1. 

A similar situation is obtained in the case of condition (18), which corresponds to 
a discontinuity in the derivative F’(x) and to a sharp angle in the curve F(z). 

4. Example. Bessel functions.” The differential equation for Bessel functions of 
order 7 can be written (L.B.H. 1, Eq. 33) 

d'f 


—; + (e* — n)f == (), (28) 
dx 


e* =z, a = logz, 
where z is the original variable. This is an equation of type (1) with 


F=e"*-n Go = +(e” — n’)™”. (29) 


’ 


The B.W.K. method can be used for 


i z2<n, x K log n n> 1, or 
(30) 
II. n KZ, x > log n 
and fails in the neighborhood of z = n, which represents a “turning point’’ (16). 
In the first interval we expand 
ee 
Gy = +in( 1 =—— 3¢ + ++) (31) 
2n 
The expansion (9) for G reads, with (10), 
G=G@,t6,+°:: =GQt = Gs Gi — 1 Ga°G3! aad (32) 


‘See J. C. Slater and N. H. Frank, Introduction to theoretical physics, McGraw Hill, New York, 
1933, pp. 148 and 347. 
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In our present case (31) it is easy to check that all correction terms can be ignored, 
since they are of order n~’, n°, ete. 


Next we compute 


a . , iz 
S = | Gdx = C + inx — e* = C + in logz — — 
Je 4n ed 4n’ 
-¢\/2 2 
om a( a+.) 
n 4n 
and our general expression (3) now yields 
\1/2 2 2 
x(=") (1 <4 254 (33) 
u= Al—] 2 ath is —¥ see ‘ 
n 4n 4n 
where K(—i/n)'’” represents just an arbitrary constant. This solution (33) obviously 


corresponds to the J, function 


“7% —_— z — eee ts « _— z — 1 eee 
J, = Ke (1 1(n + 1) T = as (1 Z(1 1) + ), 


. 1 34 
K’ = =. ( ) 
2"n! 
30th expansions (33), (34) coincide when n > 1. We thus prove that in the first interval 
(30) the B.W.IX. method leads directly to the J, solution. 
The situation is different in the second interval, where z is larger than n: 


y Qz 2\1/2 r n -2 
Go = (e" —n y' oe = aes 
(35) 
oF z n- _—s : ae no -s Y 
Go = ¢ + 9° : Gy = 4 9 = Go. 
Correction terms in the G formula (32) can no longer be neglected and 
n 3 ] l 1 
G¢= 4 a " i 
ete(-5+3 i) + i ( 4)’ 
(36) 


= lg l 2 l 
S= | Gdx = Ce +re (n’ _ 1) = c+eth(u _ ), 


u=GQ’e = a fi + (n? _ NS + - + ++), (37) 
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This reminds us immediately of the famous Debye* approximation for Hankel’s func- 
tions when z is large, or of the rather similar expressions given by Hankel: 


ee F ] es ; z 
H,’(2) = exp {-i(: _ (2n + x) (2) [P,.(z) — 2Q,(z)], (38) 
: (4n° — 1)(4n* — 9) 4n? — | 
_" a SES 7 sect = — ee 
f <7 1 2(82)" b] Q,, Rz . 


Except for a different K factor, Eqs. (37) and (38) check up to terms in 1/2’, but we 
did not push our approximation far enough in Eqs. (36) and (37) to be sure to obtain 
all terms in 1/z’ correctly. The preceding discussion exemplifies the problems discussed 
in Sees. 2 and 3 about the fact that the B.W.K. solutions obtained in different regions 
do not represent an analytic continuation of each other, and that their junction across 
the borders of these regions needs special care. Debye also developed formulas for the 
junction of his asymptotic solutions across a zero of the F(x) function; this intermediate 
solution uses Bessel functions of order 1/3 just as the Kramer solution discussed in 
Sec. 3. 
This comparison teaches us furthermore: 

1) that the series obtained by the B.W.K. method must be only semi-convergent in 
general; 

2) that they should represent a very good approximation to the actual function, since 
this is the case with the Debye series. There are, however, some fundamental limitations 
to the method of successive approximations developed in Sec. 2, Eq. (10), and a very 
striking example will be found in the discussion of the Mathieu-Hill equation (Secs. 


7 and 8). 

In such cases, it is necessary to have recourse to a different set of successive ap- 
proximations which will now be presented. 

5. A second method of successive approximations. We want to discuss Eq. (1) with a 
given function F(x), assuming we have been able to obtain the solution of a similar 


equation 
a + H(x)u = 0 (39) 


with a function H(z) that does not differ very much from the original F(z).° For in- 
stance, we may use a solution (3) of Eq. (4), where H represents the quantity in brackets 
of Eq. (4), and the G function is one of the steps of our former approximations (9), 
(10). We may try to use wu as a first approximation and write 
f=uty, (40) 
substituting in (1) and using (39) 
S’+Ff=9' + Fo + (F — Au=0. (41) 
‘Jahnke-Emde, Tables of functions, Dover, New York, 1943, pp. 137-139; Smithsonian Mathematical 
Formulae, Smiths. Inst., Washington, 1939, pp. 197-198. 


‘See S. A. Schelkunoff, Q. Appl. Math. 3, 348 (1946); M. C. Gray and 8S. A. Schelkunoff, Bell 


System Tech. J. 27, 350 (1948). 
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We now expand 


g=AnA+t+Gtgst-:::, 
gi’ + Hg, = (H — F)u, (42) 


g2’ + Hg, = (H — F)g, . 


The physical meaning of the procedure can be explained in the following way: we start 
from the wave equation (39), which obtains a solution (3) exhibiting pure propagation 
and no reflections whatsoever. We use this equation (39) to build a set of equations 
with right-hand terms representing a continuous distribution of sources along the line, 
hence an emission of secondary wavelets, propagating in both directions and starting 
from all points along the line, especially from points where (H — F) is not negligibly 
small. 

We are thus able to obtain a solution taking into account all the reflections on 
possible irregularities of the line; this should lead to a better approximation than our 
first method of See. 2. Schelkunoff’ has discussed practical methods for the integration 


of Eqs. (42). These equations read 


a d 
dae + Hg, = = h,(2x), (43) 


i = 2 (2) =H — Pgs. 
dx 


We know, from the discussion of Sec. 2, of two independent solutions u, v of the homo- 
geneous equation (39), with the condition (13). With the help of u and »v, the solution 
of (43) becomes 
ule) yf” v(x) f”’"* 
g(%:%) = —— | h,(y)v'(y) dy -— —— | h,(y)u’(y) dy (44.a) 
c y=Zo C Jy=z0 
assuming g and g’ to be zero at 2) . We thus have a method for working out the suc- 
cessive approximations (42). The h, functions represent the fictitious sources distributed 
along the line, resulting in additional waves propagating to the right (uw function) or 
to the left (v function) as explained before. 
The solution (44) can be written in a different way after an integration by parts, 


“Df nduty) dy, (44.b) 


eC de. 


- gnlerte) = — “2 [ neeyyory) ay + 
and one should not forget that g, may always contain additional arbitrary terms (au + bv) 
satisfying the homogeneous equation. 

6. The Mathieu-Hill equation, Floquet’s theorem. The preceding discussion will 
enable us to investigate successfully some types of Hill equations that were left out of 
consideration in our previous paper (L.B.H. 1). But first of all, we must restate carefully 
Floquet’s theorem; it was not stated quite accurately in the preceding paper. 

Hill’s equation corresponds to our Eq. (1) when the function F(x) is periodic. Let 
us call d the period (instead of 7 in L.B.H. 1): 


F(a + d) = F(a). (45) 
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Let u(x), v(x) be two independent solutions of Eq. (1). We obtain another set of inde- 
pendent solutions by considering u(x + d) and v(x + d), since conditions are identically 
the same at x and x + d. This means that we must have linear relations 


u(x + d) = a,u(x) + b,v(z), 


(46) 
v(x + d) = a,u(x) + b,v(z), 


expressing the new set in terms of the old one. Here we remember Eq. (13), which 


specifies that 
u(x + d)u'(x + d) — ox + djw’(x + d) = u@)o’(x) — v(x)u'(x) = C; 
hence the relation 
a,b. — a,b, = 1. (47) 


The (a, b) matrix has a determinant 1. A suitable linear combination of the original 
u, v may simplify the relations (45) by making this matrix diagonal. We choose 


U(z) = uth (48.a) 
and want to obtain 
U(x + d) = £U(z). (48.b) 
This means that 


ay, a anh = g, 


b, + b.h = th 
and eliminating h, that 
a4—-& Qa 
= —(a,+ be +1 =0. (49.a) 
b, b, —é 
The product of the two roots of this equation is 1, and we have 
t+¢° =a,+)d,. (49.b) 
The relation (48) can be expressed differently as 
U,(x + d) = tU,(z), U, = e@,(z) (50) 


t = e“*, &, periodical, period d 
and similarly 
U.(x + d) = €'U,(2), U, = e “@,(z) 


with two different functions #, , 6, both of period d. This is the general statement of 


Floquet’s theorem. 
When F(z) is an even function, it is easy to see by symmetry reasons that 


@,(x) = %,(—2z) (51) 
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in which case (and only in this case) we obtain the conditions (2) or (4) of L.B.H. 1.° 
As for Eq. (49.b), it now reads 
2 cosh (ud) = a, + b., (52) 
a general result of great importance, as will be seen in the next sections. The condition 
(47) concerning the (a, b) matrix and the final relation (52) are exactly similar to those 
obtained in the theory of filters. 

7. Stopping and passing bands, general discussion. The first essential step in the dis- 
cussion of Hill’s equation is to obtain the value of the exponent yu. A method yielding 
the uw value was discussed in L.B.H. 1 and resulted in a relation (L.B.H. 1, Eq. 18) of 
the following structure: 

cosh (ud) = y(d), (53) 
where ¥ was a certain function of d. The method developed by Whittaker results in a 
formula of similar type (L.B.H. 1, Eq. 23) and our Eq. (52) is also similarly built. The 
solution of (53) will give 

ud = Bo + ta (54) 
where £ is the attenuation constant per period d, and a, the propagation constant, 
determined modulus 27. From the form of Eq. (53) it is obvious that if a» be a solution, 
then 
al = a + 2nr (n integral) 

is another one. The function y may depend upon a variable coefficient k. In the case 
of Eq. (2), k may be the variable frequency w, or it may represent the energy £ in a 
problem of wave mechanics. 

In such problems the discussion of Eq. (53) proceeds as follows. We first plot y as 
a function of k, and obtain a curve such as the one represented in Fig. 4 (this curve 


Ad 
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a % FIG. 4 
*Since Floquet’s theorem was not used through the rest of the paper L.B.H. 1, this omission was of 


no practical consequence. 
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was drawn so as to show some typical features of the problem). We obtain passing 
bands when 
Bo = 0, or —-l<y< +1. (55) 


When, on the contrary, we obtain 
ly|>1 (56) 
this means that 8) # 0 and there is a stopping band: 
y>1, a =0 or + 2nz, 
y < -l, a = it or + (2n + 1)z. 


The a, 8 curves corresponding to Fig. 4 are plotted in Fig. 5 as an illustration of the 
general procedure. It is interesting to notice that stopping bands always correspond to 


2 = +mr (m integral) (57) 
or, if we introduce a new quantity a, 


m m 


a=-z—, a, = =-——, 
2r1 2rd 


a, being the real part of a. This condition corresponds to the general results about the 
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limits of the Brillouin zones in one-dimensional problems.’ There is, however, a funda- 
mental difference between the method followed in the present paper and the discussion 
given in the book just quoted. The method previously used was designed to give in- 
formation only about the passing bands. Propagation was assumed from the beginning 
and the » value was chosen as the independent variable 


wu = 2ria, (58) 


where a = 1/) is real, \ being the wave length. Equations (1) and (45) were then dis- 
cussed in order to obtain the corresponding value of the parameter k (k being the fre- 
quency in case of Eq. (2)). This method gave us k(a) in the passing bands and indicated 
the stopping band only indirectly by the fact that some k intervals could not be obtained 
in the discussion. 

Here, on the contrary, we start with any given k value and compute the corre- 
sponding complex uw (54) that shows either pure propagation or attenuation. Instead of 
computing k(u) we investigate the behavior of u(k), and we obtain information about 
the properties of the stopping bands in addition to the characteristics of the passing 
bands. Thus the graph of Fig. 5 shows the variation of the attenuation 8 in addition 
to the propagation coefficient a. Former graphs contained only a. 

8. Stopping and passing bands with the B.W.K. method. The important question to 
be discussed next is: how much information about passing and stopping bands can we 
obtain through the B.W.K. method? The answer is rather surprising: no practical answer 
is obtained through the direct B.W.K. procedure of Sec. 2, and only the second method of 
See. 5 can yield a reasonable solution of the problem. We gave a hint at these essential 
limitations of the original B.W.K. procedure at the end of Sec. 4. Let us explain it in so 
many words. 

The B.W.K. procedure of Secs. 1 and 2 makes use of waves (Eq. 3) that may be more 
or less distorted but always propagate freely without exhibiting any reflection. If, on 
the other hand, we investigate the cause for the stopping bands of the preceding section, 
we discover that they are due to a phenomenon of cumulative reflection: each cell of 
length d reflects back a certain amount of the incident wave, and all these elementary 
wavelets, propagating backward, happen to be sufficiently “in phase” to produce a 
large reflected wave that takes up most of the energy of the incident wave. Such a 
result can be obtained only through the procedure of Sec. 5, which takes into account 
the possible local reflections; the method of Sec. 2 cannot give any information on such 
problems. 

The discussion will be clearer if we start with a problem involving no process of 
approximation. Let us take the case of Eq. (4), the rigorous solution of which is given 
by (3). The function G(x) is any arbitrary real periodic function of x (period d). Here 
we are dealing with an equation of Hill’s type, where F(x) is represented by Eq. (5). 


The two independent solutions are 


u(a)| ioe "E iali #iS S pa [ G lr: - 
= e oOo => 7 ax; (59) 


u*(z)) 


7L. Brillouin, Wave propagation in periodic structures, McGraw Hill, New York, 1946 (p. 6 Fig. 2.4, 
p. 15 Fig. 3.9, p. 57 Fig. 15.4, p. 68 Fig. 17.2, p. 103 Fig. 27.1, p. 113, pp. 118-120, pp. 


p. 8 Fig. 2.7, 
143-145). 
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hence 
u(x + d) am Ga)" 


and 


ar+d d 
S(a + d) = S(x) + | G dx = S(x) + [ G dz. (60) 


2 


Comparing (60) and (50) we obtain 


e““ = exp (¢ | G ae) 


or 
ad 
ud = —1 | G dx + 2nzi, (n integral) (61) 
= 9 
w= —1G + ni =, 
d 
(62) 
=e 1 f? 
where G = 7 / G dx 


is simply the average value of G. Hence the B.W.K. method yields directly the very 
important quantity u, together with the fact that » is defined only modulus 7 27/d. 

The function G(x) may depend in any arbitrary way upon a parameter k, but so 
long as G remains real the formula (62) always retains its validity. Instead of the types 
of curves shown in Fig. 5 we obtain a graph of the type of Fig. 6, with curves inter- 
secting each other when a = mr. 

This proves that our Eq. (4) exhibits very exceptional features and represents a 
poor approximation to the general Hill equation. It has a continuous passing band 
and no stopping bands whatsoever. 

9. The Mathieu-Hill equation discussed with the method of Sec. 5. The discussion of 
Sec. 8 shows that the B.W.K. method does not lead to a practical solution of equations 
of the Mathieu-Hill type. The method developed in Sec. 5 will give us the necessary 
correction and enable us to obtain a much better approximate solution. 

We start with the solutions u and v = u* of Eq. (59), that satisfy a B.W.K. differ- 
ential equation (39) 

u"’ + H(z)u = 0, 
where (39) 


2 3/G’\’ | 1G” 
n=c -3(Z) +36 


according to Eq. (4). The G(x) function is periodic and exhibits the same period d as 


the original Hill equation 





f" + F(a)f = 0. (1) 


We assume that the function G has been obtained through the procedure of Sec. 2 
(Eqs. 9, 10) and that H(z) differs very little from F(x). We now use the second approxi- 
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mation represented by Eq. (44b), and we obtain two independent (approximate) solu- 
tions of (1): 
U(x) = u(x) + gi(x) = A(x)u + B(az)u*, 


(63) 
U*(x) = u*(x) + gi (x) = B*(x)u + A*(a)u*, 
with 
A(x) = 1 —- é [ (H — F)uu* dy 
; C Synz, = 
(64) 
igi 2 , 
Ba) =+[ -Pridy, = %, 
Since Eq. (43) yields 
hi(y) = (H — F)u H — F real, (65) 
for this second approximation. Our solutions (63) will satisfy relations (46) 
U(x + d) = a,U(xz) + 6,U* (2x), 
(66) 
U*(x + d) = a,U(x) + b,U*(z), 
where obviously 
a = b¥, b. = ay, 


and Floquet’s exponent results from Eq. (52) 
cosh (ud) = 1/2(a, + b.) = Re (a), (67) 
®Re meaning “‘real part of.”” Hence the only thing we need is the coefficient a; in Eq. (66). 
We can easily obtain this information from Eq. (63), taking 
zr=2+d, 
U(ao + d) = A(x + d)u(ao + d) + B(xo + d)u*(ao + a). (68) 
The discussion of Sec. (8) shows that 
u(xo + d) = exp (uod) u(2o), 
Wo = —iG 
but, at 2) we have A = 1, B = 0, hence 
u(t») = U(ao) 
and Eq. (68) reads 
U (x29 + d) = A(xp + d) exp (uod)U (ao) + B(x. + d) exp (—pwod)U*(ao). (69) 


Comparing Eqs. (66) and (69) we obtain 


a, = A(x + d)exp(ud) = exp(—iGi| 1 — 5; [- (H — F)uu* ay | (70) 
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according to (64). This gives us the coefficient a, we need in Eq. (67): 
cosh (ud) = Re exp(—iGi)| 1 v= 1] (71) 


with real 


pTot 


d 
J = | (H — F)uu* dy; 


dry 
Zo 


hence 
cosh (ud) = cos (Gd) + 1/2 J sin (Gd), (72 


and we finally obtain the Floquet exponent we were looking for. Passing bands corre- 
spond to conditions when 


—1 < cosh (ud) < +1, 
and stopping bands appear when 
cosh (ud) | > 1 (73) 


as explained in Sec. 7. The B.W.K. method of Sec. 8 never gave any stopping band. Our 
second approximation (72) may give stopping bands when 


Gd = mr + ¢, 


cosh (ud) = (—1)"(1+ 1/2 Je) +-:-. (74) 


These stopping bands appear in the neighborhood of the points Gd = mm where the 
curves of Fig. 6 intersect each other. 

This proves that the B.W.K. method must be completed with the approximation 
process of Sec. 5 in order to yield a solution of Hill’s equation. The procedure sketched 
in this final section seems to represent one of the most practical approaches to the 
problem of Hill’s equation, in addition to those cases that were discussed in a previous 
paper (L.B.H. 1). 

Thus far, we have discussed only the first step in the series of successive approxima- 
tions according to the procedure of Sec. 5. It may be found necessary to go some steps 
further and to use the expansion (42), in which case our formulas (63) could still be 
maintained. But Eqs. (64) would read 


A(z) =1- 5; ‘Me F)\(u + gi + +++ + g,)u* dy (75) 


if the n first approximations of (41) were used. This would simply mean a similar cor- 
rection in Eqs. (70) and.(71). Then 


az 


J=[ @-Putoate: + gut dy=J,+i,, (76) 
J would no longer be real and Eq. (72) should be modified to 


cosh (ud) = (1 — i) cos (Gd) + ; J, sin (Gd). (77) 


This should give a more accurate definition of the passing bands and stopping bands. 
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10. Conclusions. The principle of the B.W.K. method, as stated in Sec. 2, consists 
essentially in a procedure of successive approximations starting from a solution of 
Eq. (6), that simply represents the Hamilton-Jacobi equation of the problem: 


5So 5S5 7 7 =, 
“6x” (28 ) =f. (78) 





Go = 


In a three-dimensional problem, the wave equation (1) reads 




















Vi + Ff = 9, (79) 
and Eq. (6) yields the Hamilton-Jacobi equation 
$82)’ (25e)" (382)" - 
( da + dy * bz / . (80) 
as was shown in the paper in J. de Phys. 7, 337 (1926), Eq. 13. 
G (k) 
- \Z / 
2TT 
Tr 
) k 
if . 
-2TT /\ 





FIG. 6 


The successive approximations of the B.W.K. procedure represent a semi-convergent 
series, similar to the well-known Debye series for Bessel functions (Sec. 4). This proves 


that there is a limitation to the accuracy of the method, and that nothing can be gained 
by increasing too much the number of terms in the series. 
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The reason for this is easy to understand and was explained in Sec. 3. The Hamilton- 
Jacobi equation always leads to pure waves, propagating without reflection, as ex- 
emplified in Eq. (3). Taking the example of wave-mechanics, where 


F = 2m(E — V) (h = 2m) (81) 


we obtain in (80) the Hamilton-Jacobi equation of classical mechanics. Let us consider 
a problem of particles hitting a potential hill: they will go over the hill if their energy 
is large enough, or be completely reflected when their energy is too low, but we shall 
never obtain partial reflection with some particles reflected and some others climbing 
above the hill. This sort of thing, however, happens with the wave equation (79), and 
we discussed in Sec. 3 some typical instances where partial reflection must take place. 

This proves that after obtaining the best possible approximation with the B.W.K. 
method of Sec. 2, we must turn back to the original wave equation and use the method 
of Sec. 5. The discussion of Hill’s equation (Secs. 6-9) exemplifies the need for this cor- 
rection and shows how to apply it to a practical problem. 











381 


STUDIES ON TWO-DIMENSIONAL TRANSONIC FLOWS OF 
COMPRESSIBLE FLUID.—PART I* 


BY 
8. TOMOTIKA anp K. TAMADA 
University ef Kyoto, Japan 


1. Introduction. It is well-known that when the flow is everywhere subsonic in a field 
of flow, the nature of the two-dimensional isentropic flow of a compressible perfect 
fluid differs only slightly from that of the corresponding flow of an incompressible perfect 
fluid. Thus, in such a case, we can calculate the field of flow by any of the well-known 
methods of approximation. On the other hand, if the flow is supersonic throughout the 
field, we can determine the flow pattern by the method of characteristics. 

Difficulty occurs, however, in the calculation of transonic flows which contain both 
subsonic and supersonic regions. The difficulty arises from the fact that the fundamental 
non-linear differential equation governing the field of flow changes from the elliptic type 
in the subsonic region to the hyperbolic type in the supersonic region. 

Up to the present time, several writers have succeeded in finding exact solutions of 
the fundamental differential equation in certain particular cases, and they have thus 
obtained exact flow patterns of the compressible perfect fluid which afford valuable 
information on the complicated character of transonic flow. However, no exact solutions 
have yet been found which may be used to represent the accelerated-decelerated flow 
of gas through a nozzle. 

In the present paper, an attempt is made to clarify some of the important features 
of transonic flow. To do this, we introduce a hypothetical gas which approximates the 
real gas, subject to the adiabatic law, in the vicinity of the critical state where the 
fluid velocity becomes just equal to the velocity of sound. Then, starting from the 
equations of motion (in two dimensions) of a compressible perfect fluid, we derive a 
partial differential equation for determining, in an exact manner, the flow of this hypo- 
thetical gas. This differential equation assumes a rather simple form, but retains the 
non-linearity of the original equations. Furthermore, it changes from the elliptic to the 
hyperbolic type when the fluid velocity becomes greater than the speed of sound. 

We find that our fundamental non-linear equation can be solved exactly in several 
cases. Some of these solutions give us the flow through nozzles of various shapes. Here, 
we thoroughly discuss two flow patterns. By close examination of these flow patterns, 
we find some important and interesting features which are expected to be common to 
other similar transonic flows. 

In a subsequent paper, an alternative method of treatment will be developed on 
the basis of the fundamental equations in the hodograph plane. Just as in the case of 
the adiabatic flow of the real gas, the fundamental differential equation for determining 
the flow of our hypothetical gas becomes linear in the hodograph plane. We have dis- 
covered a group of interesting particular solutions of this linear differential equation, 
which have points of singularity of various orders in the subsonic region of the hodo- 
graph plane. Thus, by suitable linear combinations of these particular solutions we can 
obtain various flow patterns with a uniform flow at infinity. 


*Received Jan. 21, 1949. The investigations described in Parts I and II of the present paper were 


accomplished in 1946, while the main part of Part III was worked out before the end of March, 1947. 
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Although the hypothetical gas employed here has the advantage that the funda- 
mental equations are simple and can be solved exactly, it can approximate the real gas 
only for a limited range of speeds. Thus, this method has the disadvantage that it is 
applicable only to nearly uniform transonic flows. In a subsequent paper, a theory is 
developed which can be applied to flows containing limited supersonic regions in the 
vicinity of an obstacle, as well as stagnation points. For this purpose, we introduce a 
second hypothetical gas and then discuss the subsonic uniform flow past an obstacle, 
making special reference to the so-called Taylor problem: Is there any theoretical 
possibility of a continuous potential flow past an obstacle such that it flows uniformly 
at a great distance from the body and at the same time contains limited supersonic 
regions in the neighborhood of the obstacle? 

2. The fundamental equations. As is well-known, the two-dimensional steady irrota- 
tional flow of a compressible perfect fluid is given by the following simultaneous differ- 
ential equations: 


¢; = U= gq cos 6, GY =v= ¢d sin 6, 


yy. = —pqsin @, v, = pq cos 8, 


where (x, 7) are the Cartesian coordinates in the plane of fluid motion, (u, v) are the 
rectangular components of the velocity at any point in the field of flow, the angle @ 
denotes the direction of the velocity vector, and ¢ and y are the velocity potential and 
the stream function, respectively. The velocity magnitude g and the density p are 
normalized to equal unity for the critical velocity. Considering x and y as functions of 


yg and y, we then have 


l Bs 
Ly = cos 6, ty = — —sm 6, 
q pq 
a 
Yo = —sin 8, Yy = — Cos 8. 
q pq 


On eliminating both x and y, the equations of motion in the (¢, y)-plane become 


p 


1, = a(3 = La 
: p \e q 


Here, we have replaced dp/dq, using Bernoulli’s equation in differential form, namely, 


(2.1) 


dp _ _ pq (2.2) 
dq Cc b] \ 


where c is the speed of sound normalized to unity for the critical state. 

3. A theory of nearly-uniform transonic flow. In the following we shall develop a 
theory of nearly-uniform transonic flow by considering a hypothetical gas which ap- 
proximates the adiabatic behavior of the real gas in the vicinity of the critical state 
(where the fluid velocity becomes just equal to the speed of sound). 





1m Liat con a eter lineata 
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Now, let us introduce a new variable w, defined by 


“p 
w= [| 2dq. 3.1 
Le” (3.1) 


Then the fundamental equations (2.1) can be put in the form 


0, = Wy, , 
(3.2) 
6, = —Xw,, 
where 
2 
X = f, (1, a 4). 
p \q c 
If the usual isentropic law is assumed for the change of state of the gas, we have 
v++1 y¥— 1 {\i re z 
pa (plasty 
(3.3) 
gitt) tei a0 
ee ( 2 a 


where y is the ratio of specific heats of the gas. By substituting Eqs. (3.3) into the ex- 
pression for X and making use of Eq. (3.1), we can expand X(w) in powers of w: 


ie dX 44) ; 2 
= X(0) + (2 (a2 |. w+ O(w’) 


= —(y + lbw + Ow’). (3.4) 


X(w 


In the case of nearly-uniform transonic flow, the value of w is small. Consequently 
the term O(w’) on the right-hand side of Eq. (3.4) may be neglected, and the fundamental 


r«€ 


equations (3.2) can be written as 


(3.5) 


Elimination of @ from these equations yields immediately the partial differential equation 
for determining w, in the form 


(kw)yy = { (kw)*} ee ° (3.6) 


Equations (3.5) and (3.6) can only approximately represent the nearly-uniform 
transonic flow of the real gas obeying the adiabatic law. But they may be taken as the 
exact fundamental equations of motion for the flow of a hypothetical gas which obeys 
an appropriate law of change of state. This method of approach to the problem is prefer- 
able, because it is possible to discuss the field of flow with mathematical rigor. 

Before proceeding further, we shall now investigate what sort of gas corresponds to 
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the fundamental equations (3.5) or (3.6). Comparing Eqs. (3.5) with Eqs. (3.2), it is 
readily found that the relation which specifies the hypothetical gas is given by 


£23 


9 
~ 


X = —2kw, k 


Thus, taking Eqs. (3.1) and (2.2) into account, we obtain the characteristic equation 


for our hypothetical gas in the form 


{a ( L)'\ = 24°. (3.7) 
lo \pq q 























If this equation be solved subject to the conditions (p),.; = 1 and (p’),-, = —1, we 
pf C 
q c w 
q 
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Fic. 1. Solid-line curves show the curves of w, c/q and p for the hypothetical gas used in the present 
9 €/q . £ 
paper, and dotted-line curves show the curves of the corresponding quantities for the real gas 
subject to the adiabatic law. 
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obtain the curve of p(q) for our hypothetical gas, which coincides to the order of curva- 
ture at the critical state gq = 1 with the corresponding curve of p(q) for the real gas 
subject to the adiabatic law. The curves of p(q), c(q), and w(q) for our hypothetical gas 
when y = 1.4 are shown by solid-line curves in Fig. 1. The corresponding curves for 
the real gas obeying the adiabatic law are shown by dotted-line curves, for the sake of 
comparison. 

4. Some exact solutions of the fundamental equation (3.6). The fundamental equa- 
tion (3.6) for our hypothetical gas can be solved exactly by assuming kw in the following 


forms: 


(a) kw = fie + Ap), (A = a certain constant) 
(b) kw = &(¢) + W(y), 

(c) kw = B(¢g)VW(y), 

(d) kw = W(y) + W(We’, 


Zot v) + 2y, 


where Z(¢ + yw’) denotes a function of g + y’. 

These solutions together with the corresponding flow patterns have been investigated 
in detail by the junior writer in his doctoral thesis.’ Only a brief sketch of the results 
obtained for the first three cases will be given here, however. The remaining two cases 
will be discussed in some detail in the following sections. 

The first exact solution (a) gives the so-called spiral flow. Since an exact solution 
giving the similar flow pattern for the fundamental equations of motion of the real gas 
is well-known, we can compare the two solutions. It has been found that the flow of 
our hypothetical gas approximates sufficiently accurately the flow of the real gas. 

The second solution (b) represents the flow through a Laval nozzle (i.e., a converging 
and diverging nozzle), which changes from the subsonic state to the supersonic one as 
it passes through the narrowest cross section. The results have been compared with 
those of the approximate calculation due to Meyer’ for the corresponding adiabatic 
flow, and excellent agreement has been found. 

The third solution (c) also gives the flow through a nozzle having a straight axis. 
Detailed calculations have been carried out for the flow which contains both the curves 
of singularity, / = 0 and J =o, where 


A(x, y) _ 0(x, y) Ay, y) z 


~ aq, 8) ay, W) aq, 8) 


5. The flow through a nozzle with two contractions. Let us consider the fourth solu- 
tion (d), of the form 


™ 
| 





kw = V(y) + Uwe’. (5.1) 


1K. Tamada, Studies on the two-dimensional flow of gas, with special reference to the flow through various 
nozzles, to be published. 
Th. Meyer, Uber zweidimensionale Bewegungsvorgdnge in einem Gas, das mit Uberschallgeschwindig- 
keit strémt, Dissertation Goettingen, 1908. 
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Substituting this into the fundamental equation (3.6), we obtain two differential equa- 
tions for determining the functions ¥(y) and W,(y). After considerable manipulation*, 
the solution for V, can be shown to be 

VY, = Ae(t + w.), ¢ = (2r)'"y, 


where A and w, are constants of integration, and ¢ is the,Weierstrass elliptic function.* 
The quantity w, is chosen to equal the real half-period of,the g-function so that the 
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Fic. 2. Case when g3 = —1. 


























flow may be symmetrical about the central streamline y = 0. With this result for Y, , 


the function VY, becomes 
ail ae" ( = dey" 
== ai(1 Tar tS 


lg 


or 
Vv, = 2°°ao'” cos (3 arctan — ?) 
. S 7 dt 


according as the invariant g; of the g-function is taken equal to —1 or +1. 





*Since this flow will be of less technical interest than that of case (e), we omit the manipulative 


details. 
3For the Weierstrass elliptic function g used here, see E. Jahnke and F, Emde, Funktiontafeln mit 


Formeln und Kurven, Teubner, Leipzig and Berlin, 2nd ed., 1933, p. 166. 
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Substituting the expression for w as given by (5.1) into the fundamental equations 
(3.5) and carrying out a simple calculation, we obtain the result that 


it son dV, ld¥, 5 (5.2) 


dy** 3 dy*’ 
where an additive constant of integration has been omitted since @ is everywhere zero 


along the streamline y = 0. 
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Fia. 3. Case when g3 = +1. 























The curves of VY, and W, are shown in Fig. 2 for the case g; = —1 and in Fig. 3 for 
the case g, = +1. From these figures, in conjunction with Eqs. (5.1) and (5.2), it will 
be seen that the present solution gives different types of flow depending on the choice 
of signs for the constants g; and a. 

Detailed numerical calculations have been carried out for the case gj = —1, by 
assuming a = 0.15, \ = 0.5 and taking, as usual, y = 1.4 for air. Figure 4 shows the 
flow pattern in the physical plane and Fig. 5, the corresponding flow pattern in the 
hodograph plane. The correspondence between the physical and the hodograph planes 
can be easily understood, if we regard the streamlines in the hodograph plane as a group 
of contour lines forming twisted surfaces. 

The flow pattern in the physical plane does not show any peculiarity, but, as is 
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Fig. 4. Physical plane. Heavy solid-line curves and thin full-line curves represent streamlines and 
equi-velocity lines, respectively, while dotted-line curves represent Mach waves. 
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Fic. 5. Hodograph plane. Heavy solid-line curves represent streamlines, while dotted-line curves 
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clearly seen from the flow pattern in the hodograph plane, it contains the particular 
Mach waves AB, A’B’,CB and C’B’ as the curves of singularity J =o. 

6. Transition from the flow of Taylor’s type to that of Meyer’s in a given Laval 
nozzle. As the last example, we shall consider the flow represented by the fifth solution 


(e), of the form 
kw = Z(u) + 2y’, (6.1) 
where 
u=ety. 
Inserting this expression for w into the fundamental equation (3.6), we obtain the 
following equation for determining the function Z(u): 


d (z a) _ dZ _ cae (6.2) 
du du du 


This equation can be integrated once immediately, and we get 
{ 5 Af) 5 


Z a tn teak (6.3) 


where a constant of integration has been included in wu for the sake of convenience. 
If we set dZ/du = xX; it follows from Kq. (6.3) that 


u= ; (x — 1)Z. (6.4) 


=5 


Further, again taking into account that dZ/du = x and d’Z/du® = xdx/dZ, we have 
= xX ’ 


from Eq. (6.2) that 


r dx 9 
¥ —— —9 ). 
Zx dZ (x (x 1) 


which can be integrated immediately to give 
Z = a(x — 2)°°(x + I)”, (6.5) 
where 2°°a is a constant of integration. Eliminating x from Eqs. (6.4) and (6.5), we get 
ultimately 
(Z — 2u)?(Z.+ u) = 2a’. (6.6) 


Next, inserting the expression (6.1) for w into Eqs. (3.5) and carrying out simple 
integrations, we obtain the expression for @ in the form 


/ —— 
k= au z + 2¢+ 3 v), (6.7 


where a trivial constant of integration has been neglected. 
From Eqs. (6.1), (6.6), and (6.7) it is seen that the flow under consideration is 
symmetrical about the straight streamline y = 0. Therefore, replacing two streamlines 


) = —y, by two solid walls, we obtain a Laval nozzle having the line y = 0 


yp y, and y 
as its axis. 


Now, it is well-known that there may be two different types of flow in a given Laval 
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nozzle. In one type, the flow is subsonic on both sides of the throat of the nozzle but 
contains limited supersonic regions in the neighborhood of the walls at the throat 
(Fig. 6), while in the other the flow changes from a subsonic to a supersonic state as it 
passes through the throat (Fig. 7). The former type of flow was first discussed by Taylor 
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Fic. 6. Flow of Taylor’s type. Fia. 7. Flow of Meyer’s type. 


and therefore will be called ‘“Taylor’s type’’. On the other hand, the latter type of flow 
was first investigated by Meyer in a paper already mentioned, and hence it may be 
conveniently referred to as ‘‘Meyer’s type”’. 

It is expected that under certain circumstances the transition from Taylor’s to 
Meyer’s type may occur in the flow through a Laval nozzle. The problem has already 
been investigated to some extent by Taylor‘, Gértler’ and others. Unfortunately, how- 
ever, these authors have not been able to arrive at any definite conclusion, since they 
have employed only approximate solutions in the form of series in x and y, which be- 
come rapidly divergent when the transition state is approached. It will be seen subse- 
quently that it is possible by our present solutions to discuss the problem of transition 
from the flow of Taylor’s type to that of Meyer’s in a given Laval nozzle. 

First, the values of the function Z(u) must be calculated by using Eq. (6.6) in the 
vases a > 0 anda < 0. The curves of Z(u) plotted against u are shown in Fig. 8. These 
curves represent the velocity distribution along the central axis y = 0 of the nozzle, 
since the function w, which is connected with g as shown in Fig. 1, becomes equal to 
Z when y = 0. It can be seen from these figures that the case a < 0 is suitable for our 
the problem of transition of flow from Taylor’s to Meyer’s 


purpose of investigating 
In the case a < 0, the curve of Z consists of two branches A 


type in the Laval nozzle. 
and B as indicated in Fig. 8. On the branch A, Z is always negative and becomes equal 


to —~ as u —+-. Therefore, kw = Z + 2y’ is everywhere negative (subsonic) for 
small values of | y |, but it becomes partly positive (supersonic) in some limited portions 
on those streamlines for which the values of 2y° are larger than | Z,,.. | = | a|. Thus, 


the solution in this case would give a flow pattern of Taylor’s type as shown in Fig. 6. 
As for the branch B, we first consider the range in which Z < 0. Then, it is seen 
from Fig. 8 that both Z and wu increase as ¢ increases along any streamline y = const., 


‘G. I. Taylor, The flow of air at high speeds past curved surfaces, British A. R. C. Reports and Memo- 


randa, No. 1381, 1931. 
5H. Gortler, Zum Ubergang von Unterschall- zu Uberschallgeschwindigkeit im Diisen, Z. angew. Math. 


Mech. 19, 325-337 (1939). 
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until the value of dZ/du becomes infinite at the point u = uw), Z = 0. This means that 
the flow has reached the curve of singularity J = 0; for we have 

O(q, 8) — d(q, 8) (yg, W) 

A(x, y) Ay, W) A(x, y) 


_ (2° (3 7 a a oe 
= 3 1" 5 oe g Ve qv f = s 


since q, and q, are proportional to w, and wy, and hence contain dZ/du by Eqs. (3.1) 


J = 


and (6.1). Therefore, the flow in this case would be as shown in Fig. 9(a). 


\ \x7t> 


. 


















(b) 


Fic. 9. 


0 of the same branch B, we obtain a flow pattern as shown in 


For the range Z > 
= (0 (the 


Fig. 9(b), which is everywhere supersonic and has the curve of singularity J 
shock line) in common with the preceding flow in Fig. 9(a). 

In the case a = 0, the branch A coincides with the two asymptotes a0 and )0. In 
this case, the local speed of sound is attained at some point (which corresponds to the 
point 0 in Fig. 8) on the axis of the nozzle; hence the upper and lower limited supersonic 
regions come in contact with each other at that point. This may be considered as the 
limiting case of the flow of Taylor’s type. The flow pattern is as shown in Fig. 10(a). 

In this case (a = Q), the branch B is also coincident with the asymptotes bO and 
cO (=a0) and hence the two branches A and B may be joined at O into one straight 
line aOc. If this branch aOc is adopted as Z, our solution would obviously give a flow 
pattern of Meyer’s type as shown in Fig. 10(b), which has part of the flow pattern in 
10(a), corresponding to the half-branch aO 


common with the preceding case in Fig. 
» transition of flow from Taylor’s to Meyer’s 


of Z. Figures 6 and 10 suggest a possible 


type in the Laval nozzle. 
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Taylor and Gértler have treated the flow between two rigid circular walls symmetrical 
with respect to the z-axis. The dimensions of their nozzle were determined to satisfy 
the condition 

R = 4h (6.8) 
at the throat of the nozzle, where R denotes the radius of curvature of the walls’ and 
h is the half-breadth of the throat. 


mis, ana a 
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Fia. 10. 


For the sake of convenience, we shall use the same condition to specify the shape of 
the nozzle.’ This condition determines the particular streamlines to be chosen as the 
walls of the nozzle, for given values of the parameter a. As will be seen presently, these 
streamlines are almost coincident over a considerably wide range of zx, for the several 
values of a. Thus, our nozzle retains approximately a fixed shape in spite of the change 
of flow, and, therefore, we can discuss the problem of the transition from Taylor’s to 
Meyer’s type of flow in a given Laval nozzle. 

Now, at the point (go , Yo) on the wall at the throat of the nozzle, @ is obviously 


zero and hence we have, from Eq. (6.7), that 


9 
Z= —-Wy — 3 Vi. (6.9) 
Remembering that y, = cos 0/pq the half-breadth of the throat of the nozzle is given 
by 
evo 1 
h=| (—cose dy. (6.10) 
"oe pq =o 


Now, the radius of curvature R of any streamline is given by 


_ds _dsdgp_ 1 
~ d0 dedé~ q(d0/dy) * 





R 


‘In Taylor’s and Géortler’s cases, R is the radius itself of each circular wall. 
7It may be remarked here that unlike Taylor’s and Gértler’s cases, the walls of our nozzle do not 


assume the exact circular shape. 
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By: Eqs. (6.7) and (6.3), we have 
d@ 2Qy (s au) 
de k\ + ZI 
Therefore, 
kZ 
TAF LOA (6.11) 
2q¥(3Z + 2u) 


2 


t= 


Thus, taking into account Eq. (6.9) the radius of curvature at the point (go , Yo) 
on the wall is given by 
k(3¢o 0) 
ee (6.12) 
12 qo¢0 Wo 


Yo). In the present case, if use is made of Kq. 


R = 


| 


where q is the value of q at the point (¢ , 


(6.10), the condition (6.8) can be expressed in the form 
k(3¢o 0) oad &. 
x + 5 = 4 | (. cos 0) dy. (6.13) 
12 gogo Wo “0 Pd P= Po 


If the two equations (6.9) and (6.13) are solved by the method of successive approxi- 
mations, the values of g) and y for various values of a can be determined. We have 


y 
h 
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calculated the values of ¢, and y, for three cases in which a = —0.10, a = —0.05 and 
a = 0, respectively, obtaining the results as shown in Table I. 


Table I 


a —0.10 —0.05 0 
Yo 0.05028 0.04665 0.04605 
Wo 0.3620 0.3707 0.3717 
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Figure 11 shows the shapes of the streamline y = y , which is to be replaced by 
the wall of the nozzle, as calculated for these three cases, by taking 2/h and y/h as 
abscissa and ordinate, respectively, and choosing the point (go , %o) as the origin, the 
ordinate being shown on an enlarged scale. From this figure it will be seen that there is 
very satisfactory coincidence between these streamlines and, moreover, that if we make 
adjustment by shifting slightly the respective origins for the cases a = —0.10 and a = 
—0.05, almost perfect coincidence can be obtained within the range of 2/h shown. 

Now, as already mentioned, when a = 0 our solution gives two flow patterns, namely, 
the limiting case of the flow of Taylor’s type and the flow of Meyer’s type, which in 
the upstream subsonic region are exactly the same. The dotted-line curve in Fig. 11 
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Fig. 12. Heavy solid-line curves and thin solid-line curves represent streamlines and equi-velocity 
lines respectively, while dotted-line curves represent Mach waves (characteristics): a = —0.10. 
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shows the shape of the corresponding streamline for the case of the flow of Meyer’s 
type. It is seen that the deviation of this curve from the group of curves corresponding 
to the case of the flow of Taylor’s type becomes rather conspicuous in the supersonic 
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region on the left-hand side of the throat. However, this fact seems to have little effect 
on the discussion of the transition of flow, because in such a supersonic region, the wall 
can be continuously deformed by a small amount without causing any appreciable 
change in the upstream subsonic region. 
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We have calculated the flow patterns in the physical plane as well as in the hodo- 
graph plane, for the three above-mentioned cases, namely, a = —0.10, a = —0.05, 
and a = 0. The results are shown in Figs. 12-15. Also, the velocity distributions along 


the axis of the nozzle are shown in Fig. 16. 
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Comparing these figures we can explain the possible transition of flow from Taylor’s 
to Meyer’s type in a given Laval nozzle in the following way. As the rate of flow in- 
creases, the curve of velocity distribution along a streamline becomes steeper until, in 
the limiting case of the flow of Taylor’s type, it reaches an infinite curvature on a partic- 
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ular Mach wave AB (the curve of singularity J = ©) as shown in Fig. 14. In this stage, 
the local speed of sound is first attained at some point A on the axis of the nozzle. At 
the same time the streamlines in the hodograph plane come in contact with a particular 
characteristic curve AB which corresponds to the Mach wave AB in the physical plane. 

In the next stage, if the conditions governing the flow (such as the pressure at the 
outlet of the nozzle) be changed sufficiently, the flow changes abruptly to that of Meyer’s 
type in the downstream region, having the Mach wave AB as the natural boundary of 


. , 
q 











— 1.1 
/ 
a=0 
Meyer's type 
1.0 














(~~ 
tg > | > ~ 0.9 
een in. 



































ZY ~~ ~a=-005 } Taylor's type 
ees 
-0.8 -0.4 A 0 04 X/h 08 
Fic. 16. 


the unchanged upstream region. Such an abrupt change seems to be inherent in the 
flow, since the condition that the streamlines be in contact with some particular char- 
acteristic curve in the hodograph plane, in conjunction with the fixed shape of the 
wall of the nozzle, seems to determine the unique and isolated solution for the flow of 
Meyer’s type. Moreover, the flow of Meyer’s type cannot be approached in a continuous 
manner from the group of solutions for the flow of Taylor’s type. 

If the change of conditions controlling the flow is not sufficient for the occurrence 
of such an abrupt change of state, then intermediate stages occur in which the flow 
forms shock waves in the downstream side so as to fit the conditions at the outlet of 
the nozzle. At the beginning of these intermediate stages a slight shock wave first 
appears along the particular Mach wave AB (the curve of singularity J = ~). 
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ON THE AXIALLY SYMMETRIC FLOW AROUND A 
NEW FAMILY OF HALF-BODIES* 


BY 
A. VAN TUYL 


Naval Ordnance Laboratory, Silver Springs, Md. 


1. Introduction and discussion. A recent paper by A. Weinstein [1]f includes a dis- 
cussion of the flow of a perfect incompressible fluid corresponding to a uniform disc- 
source at right angles to a uniform parallel flow from infinity. There is a stream surface 
in the form of a half-body under certain conditions whose shape depends on a single 
parameter. As the radius of the dise goes to zero, with both the total source strength 
and the velocity at infinity held constant, the half-body continually decreases in blunt- 
ness until the form of the well-known Blasius-Fuhrmann half-body (cf. [1], for example) 
is reached. When the radius increases under the same conditions, an upper value is 
reached at which the nose coincides with the disc. 

The object of the present paper is to make a combined numerical and analytical 
study of these half-bodies, and some features of the corresponding velocity distributions. 
A set of graphs is obtained by means of which the nose of any given profile can be 
plotted; representative profiles are shown. Curves are also given which show the position 
and the magnitude of the maximum velocity on the surface as functions of the half- 
body. It follows from Bernoulli’s theorem that the pressure is least at the point of 
maximum velocity, hence the danger of cavitation is greatest there in an actual flow. 
All calculations are made using expressions for the stream function and velocity compo- 
nents in terms of elliptic integrals. 

The total source strength will be denoted by m, the radius of the dise by b, and the 
velocity at infinity by U. Cylindrical coordinates x, p will be used with x measured in 
the direction of the parallel flow and with the center of the disc at x = 0, p = 0. Then 
as given in [1], the stream function for the total flow is 


- (0, z<0 
vy = 5 Up F 2mpb"' | e'7'*J,(pt)J,(bt)t"' dt + 4 (1.1) 
’ ™ | —2m, xr>0 


where J,(z) is the Bessel function of the first kind of order one. The additive constant 
is chosen so that the stream function vanishes on the line p = 0, x < 0. It is shown in 
[1] that b < (2m/U)"” is the necessary and sufficient condition for y = 0 to be the 
equation of a half-body. The b-independent quantity (2m/U)'” is the intercept of the 
profile on the p-axis, and will be denoted by p) . The half-body is asymptotic to the 
cylinder p = 2'"p,. 

The components of velocity in the z- and p-directions, respectively, are 


_, OW _; OW 
1 1 
“= —_ v= — <i 1.2 
Pap’ i (1.2) 
*Received March 29, 1949. 
tNumbers in brackets refer to the references at the end of the paper. 
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Concerning the stream function for the uniform disc-source, some errors have oc- 
curred in the literature. The stream function was given incorrectly by Beltrami [2], 
who apparently overlooked the many-valued character of the stream function for the 
ring and integrated a branch which is continuous across the interior of the ring. This 
error was reproduced by Webster [3, p. 367]; Lamb [4, p. 239] has given it incompletely, 
but correctly as far as it is given. The flow corresponding to a ring of sources has been 
investigated by several others. However, it seems that the fundamental fact that the 
corresponding stream function is many-valued has not been mentioned except in [1]. 

We now make the substitution bt = s, and introduce the notation 


y*=m'y, (\=b'p, E=dD'2, n= db'p, u®* =U 'u, v* = U's. (1.3) 


For (1.1) and (1.2), respectively, we obtain 


2 200 0, é<0 
y* = um F 2 | e'*!* J. (ns) J,(s)s* ds + 4 (1.4) 
r Jo ” J 
—_—s £ > 0 
’ 5 2. 
rn? dy* ” dy* 
u* =- ’ rs —- ., (1.5) 
2n On Zn o& 
From (1.4) and (1.5), it follows that 
- (fe<0 
u*=1+FH | e'*'* Jo(ns) J,(s) ds, 4 (1.6) 
° le>0 
y* = }? | e'*'*J.(ns) J,(s) ds. (1.7) 
Given \, the equation of the profile Y = 0 in the &y-plane is 
= (0, t<0 
Nn? = +277" | e~'*'* J,(ns) J,(s)s"' ds + 4 = T(E, »). (1.8) 
™ \2n", §&20 
We introduce also 
ut = 1+ (XT) "'(u* — J), vt = (XT) 'v*, we = (u®? + v**)!”? (1.9) 


3y 1.8, w* at a point (€, 7) is the velocity magnitude along the profile which passes 
through it. 

It is possible to express (1.6), (1.7), and (1.8) in terms of elliptic integrals, and this 
is done here. It was noted by Weber [5] that integrals of the form of (1.8) can be evalu- 
ated in terms of elliptic integrals, although the reduction was not carried out. Both 
complete and incomplete elliptic integrals of the first and second kinds are present in 
(1.6) and (1.8), while the evaluation of (1.7) involves only complete elliptic integrals. 
As earried out here, (1.7) and (1.8) are evaluated separately, and (1.6) is expressed 
as a linear combination of (1.7) and (1.8) and a third integral which can be evaluated 
in terms of complete elliptic integrals. We thus arrive at expressions for (1.6) and (1.7) 
of the form a,K + 6,EF + cT + d and a,K + bE, where K and E are the complete 
elliptic integrals of the first and second kinds, and the coefficients are simple algebraic 








1950] AXTALLY SYMMETRIC FLOW AROUND HALF-BODIES 401 


functions of £ and ». It turns out that (1.8) is most simply expressed in terms of K and 
(K — E)/k D, where k is the modulus, and (1.6) and (1.7) in terms of K and 
(2D — K)/k’ C. Both C and D are tabulated in Jahnke and Emde, [6]. The tables 
of K and # are the most accurate, however, and are used here throughout. 

In the first part of the computations, 7'(é, 7) is calculated as a function of & for several 
fixed values of », and the results are plotted. Then making the substitution ¢7’”* = 2, 
aaa p, a second set of graphs is made from the first showing’ 7 as a function of x 
for several fixed values of p between 0 and 2'”. It is evident that the profiles T(z, p) = 


constant all have the p-intercept unity, and are asymptotic to p = 3”. In Ar pea 
(9T'/dx), has been found as a function of p at T = 0, for 0 < p < 2'”, wl at T = 1, 
x 0, where 0 < p < 1. The first is a simple algebraic expression, the second involves 
elliptic integrals. Figure 1 shows 7' as a function of x for several values of p, while Fig. 2 
gives some of the profiles 7’ = constant. 

To investigate the position of maximum velocity along the profile, the function 
R(é, 7) 2°'u* (@w**/dE)y has been obtained in terms of elliptic integrals. It is ex- 
pressed in the form au** + bu*v* + co*¥*, where the coefficients are linear combinations 
of complete elliptic integrals and 7'(é, 7). Evidently, the curve R(é, ») = 0 intersects 
a profile 7 = constant at the point of maximum velocity. To find points of R(é, n) = 0, 


u* and v* have first been calculated as functions of — for the same values of 7 taken 
before; then the zero of R(£é, 7) has been found for ial n. The curve begins at & = 0, 
n | and crosses the axis at some 7 > 1. This value of 7 has been found numerically. 
From the graphs of 7'(é, 7) and the transformation 2 = £7’, the value of x at which 
maximum velocity occurs, denoted by z,, , is found as a function of 7’. It is found that 
%m(7') is positive forO0 < T < 0.819, and negative for 0.819 < 7 < 1. At T = 0, which 
defines the Blasius-Fuhrmann half-body, z,,(7’) is equal to 6 '”*. Further details are 
that dx,,/dT = —(13/48)6'” at T = 0, and dz,,/dT = +0 at T = 1 — 0. When meas- 
ured from the nose instead of from the cross section of radius unity, 2, seems to be 
always decreasing. The plot of x,, as a function of T is shown in Fig. 3. The portion of 


the curve between 7 = 0.819 and T' = 1 is below the axis, but too close to be distin- 
guished from it. 

Finally, w*nax iS plotted as a function of 7 in Fig. 4. It becomes infinite at 7 = 1 
and is equal to 2/3'” at T = 0.* 

2. Evaluation of 7'(£, 7). Setting «4 = 0, v = 1 in Watson [7, p. 389, formula (1)], 


we obtain 


, 


Loe Sint J (Ot dt = (2x?) [ F(1, 3/2; 2; —p’t~’) sin’ ¢ de, (2.1) 


“0 
where F(a, b; c; z) is the hypergeometric function, and p’ = 1 + n° — 27 cos y. We have 


: } soe (ie — FD as 
f ¢ ).- D+» = <n 
F(1, 3/2; 2:2) = 2 ) er (Qn) "i 


| (2.2) 
= 2°((1 — 2)? — 1). 


*The author is indebted to Miss C. Brudno and Miss G. Fulton for the drawings and numerical 


calculations. 
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Thus, 
Osi ded a sin’ yg dy 
"T(nt) J (it dt = 4 | FF _ 
I . deities : wrJo l+n — 2n cose 
(2.3) 
_ én 4 se sin® ge dg = 
tr Jo (1+ — 2n cosy) +1 + n° — 2n cosy)” 
We find that 
| (1 + n° — 2n cosy)" sin” gdp = 553, n> 1, 
(2.4) 
T o 
= 2 ’ ns l 


— 1 into the second integral on the right and set 4n[é° + 


We substitute cos ¢ 
)~? = P. Clearly k’ < I? < 1. Then for all > 0, 


= 2 
(n+ 17)" =k’, 4n(1 + 9 
[ (1 + n° — 2n cosy) '(& + 1+ 7° — 2n cosy)’ sin’ g dp 


70 


el 
7 fa- aa Rey dt (2.5) 


“0 


= 81+) +0 + ol 


— (I? — 1) | fa—-PP) "da - ’/) 7a -— Rey” at|. 


Finally, to reduce the integrals on the right to standard elliptic integrals, we substitute 


t = sn (u,k). We have* 


al K 
| 1 -— P71 — RP) dt = [ sn’ (u,k) du = (K — E)/k’, (2.6) 


“0 *0 


where K and E are the complete elliptic integrals of the first and second kinds, re- 
spectively. Similarly, 


— KR?) dt = [ (1 — I? sn? uw) sn’ u du. (2.7) 


) 


pl 

| (1 we rey = P7701 
With sn (6, k) = I/k, 

K oK 7,2 2 

> 9 _— | 9 9 = “g S rs 5 S 1 , - 

[ (1 — I’ sn’ u)7' sn? u du = (k’ sn 6 cn 6 dn 8)7’ | ken $n 6 dn 3 7 pa . (2.8) 
Jo Jo 1 — k*sn* 6sn° u 


The integral on the right of (2.8) is II(K, 6), the complete elliptic integral of the 
third kind with parameter 6, having as its path of integration the segment of the real 
axis from 0 to K. We note that II(K, 6) is doubly periodic in 6 with the periods 2K 
and 2K’i. By [8, p. 523] we have 


II(K, 6) = nxi + KZ(6, h), (2.9) 





*All formulas concerning elliptic functions and elliptic integrals used here are given in Whittaker 


and Watson [8, Ch. XXII]. 
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where Z(6, k) is the Jacobi Zeta function, and n is an integer which depends on the 
determination of 6. . 

We now make the substitution 6 = K + iM, and use the addition theorem, obtaining 
Jacobi’s imaginary transformation 





_ en (iM, k) _ 1 
a= ab" aun” (2.10) 


where k’”” = 1 — k’. Since k < k/l < 1, it follows that M is a quantity between —K’ 
and K’, plus multiples of 2K’. Noting that II(K, 6) and Z(é, k) both tend to zero as 
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° 
-0.8 -0.4 0 x 0.4 0.8 1.2 


Fic. 1. 


M tends to zero, we find that n = 0 for —K’ < M < K’. We see that II(K, 6) has 

discontinuities at M = (2m + 1)K’, and that II(K, 6) = 0 at M = (2m + 1)K’. The 

discontinuous behavior of II(K, 6) can also be seen from (2.8), since the left side be- 

comes infinite at both M = K’ — 0and M = K’ + 0, while cn 6 changes sign at M = K’. 
All values of M differing by multiples of 2K’ are given by 


M = [ (1— k? sin? 9)" de (2.11) 
Jo 


when ¢ = are sin sn (M, k’) = are tan [¢(7 — 1)~’] runs through all values differing 
by multiples of x. Choosing ¢ so that ¢ = 0 for = 0, 7 > 1, we see that —K’ < M < 
K’ when » > 1. As a function of g, M is Legendre’s incomplete elliptic integral of the 
first kind with modulus k’. In the notation of Legendre for incomplete elliptic integrals, 
used in Jahnke and Emde, 


M = F(k’, ¢). (2.12) 
By [8, p. 518], 
Z(6, k) = E(6, k) — (E/K)é, (2.13) 
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where (6, k) is the incomplete elliptic integral of the second kind with argument 6. 
Using the addition theorem and Jacobi’s imaginary transformation for #(é, k), we 
obtain, for —K’ < M < K’, 


— — J p2sn CM, k’) en (M, &’) a al le i 
IM(k, 6) = i iP + (K — E)M — KE(M, | |. (2.14) 


The notation is that of Whittaker and Watson, who write 


aM 


E(M, k’) = | dn? (u, k’) du. 


“0 
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Fic. 2. 
In the notation of Legendre, 
E(M, k’) = E(k’, ¢), (2.15) 


where M and ¢ are connected by (2.11). 
Using the addition theorem in each case, we see that cn 6 is negative imaginary, 


while dn 6 is positive. Thus, 
sn 6 = L/k, en 6 = —ik'(l? — k’)'”, dn 6 = (1 — T°)” (2.16) 
and 


k’ sn 6en 6dn 6 = —il(l — 1°)'?°(P? — k*)™. (2.17) 
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It follows from (2.3) that 
(0, e<0 
r=4 foc! SQA (noe at ae 
is continuous at £ = 0, 7 > 1 having the value (27)~*. Then for 7 > 1, substituting back 
and replacing K — E wherever it occurs by k*D, we obtain 


—&2>0 


ST 


into (2.5) and (2.3) 


oe. 3 |p 4 2-1 K | +21 Ke prRU, 0 — KE’, ¢)]. (2.18) 
4n 2r7n 


Since both sides vary continuously through 7 = 1, (2.18) must hold also for 7 < 1, 
£ ~ 0. Finally, substituting into (1.8), 


9)-& oe 2 2 _ Z 7 
+ =, E +- (1 ) K | + us — (k°DF(k’, ¢) — KE(k’, ¢)], (2.19) 
TN 4n TN 


1, 2) = 

” 
where » = are tan [£(n — 1)~'], with ¢ = 0 at € = 0, > 1. In the region outside the 
disc, ¢ runs from —7 to x. For purposes of calculation, we have the following less 


compact form of 7'(é, 7): 


7=> I] 
T(E, n) = 4 A | D+ 7 D K| 
n T } 
hat (e<o 
= 1} [k’DF(k’, @) — KE(k’, 2) < 
sa :>0 
(2.20) 
7 7 ] 
oy A — 2 
T(t, ») = 414+ 2% | p ,2=—2 K| 
™ 7) 
(0, e<0 


+ uae ‘DF (k’, @) — KE(k’ ¢)] + 4 
si \Qn*, €>0 


where @ are tan | &(7 — 1)7' *g E. o<T 2. We obtain (2.20) from (2.19) using the 
relations F(k’, 2x + ¢) = 2K’ F(k’, ¢), E(k’, 2x + v) = 2E’ + E(k’, ¢), and the 
Legendre relation KE’ + EK’ — a =f . 

AtT =0,0< p< 2™, 


aT’ is : — 
(22) = — 2 (p? — 12 - p)”. (2.21) 
Ox/ » 3 


This result is most easily obtained by going back to (1.8) and (1.12) and making the 


substitution £ = 27" '”*, » = pT 
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As seen in Fig. 1, the point 7’ = 1, x = 0 is the intersection of all curves p = constant 
for 0 < p < 1. At this point, 
$f) “1 -3 
=) « i Y 2.22 
(2) 16x C(1 + pp)”, ( ) 
where k? = 4p(1 + p)~”. 
3. Calculation of the velocity components. By partial integration, noting that 
dJ,(Z)/dz = Jo(z) — Ji(z)z7", we have 


n | e*S,()[Jo(nt) — Ji(nt)(nt)'] dt = / 1) 2 Ind dt 


=— [ ® Int -B0 + J,(i) -—Ji()t']dt, ¢€>0. 


“0 


Further, 


nf ce S(ICnt) at = — [eI & Jolnt) at, 


“0 


II 


1+ [ &*Jo(nl—eso(t) — Ju(o) at. 
Thus, we establish the identity 
Qa) [ec ®S(aI) dt =1-29 | eos or at 
(3.1) 
—en [| cad sy at —& | "Ion Jolt) at 


By Watson [7, p. 389], 


a co 


| ec J,(nt) S(t) dt = Q(t FEE), 


where Q,(z) is the Legendre function of the second kind of order v. For v + 1 > 0, we 
have the integral representation (cf. [8, p. 317]) 


Q,() = 2°77 [ (1 — fe -— 7°" dt; 


substituting ¢ = 2u? — 1, with 2(1 + z)~’ = k’, into this, 


1 
Q_1(2) =k [ (1 — wv) "7(1 — ku’)? du = kK, (3.2) 
Qin) = fw — wy — ety? du. (3.3) 


When z = (1 + & + 7)(2n)7', we have 2(1 + z)~* = 4n[é? + (1 + »)’]"’, the quantity 
already defined as k*. Substituting u = sn (v, k) into (3.3), we find that 


K 2 2 
Qinle) =F [ De ae, (3.4) 
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We may easily verify the relation 


> fsn’u en’ u sn u en u ‘ 
k [= — du = —u + —— 42 f sn u du, 


dn* u dn u 


from which it follows that 
(2) = KC. (3.5) 


By (1.8), 


si (0, é€<0 
Qn f cmd (Ot dt = 40°, ») +4 
: 2, €>0. 


Substituting the preceding results into (3.1), with £ > O replaced by | ~ |, we get 


< 0, 


wtr 


: aw sa oa ke é . | 
(l — 9) ee J o(nt)J() dt = 41 F of TE, 9) + —ia (K + FC) 4 
Jo Tn le > oO 


ctr 


Finally, substituting (3.6) and (3.5) into (1.6) and (1.7), we obtain the following ex- 
pressions for the velocity components: 


\7(u* — 1) = (1 — 7) ‘In°T(é, n) — ké(arn' | "(K + k’ nC) — fT), (3.7) 


N~*y* = h(a”) 'C. 


-] o 2/a:¢ * 
2 ‘ux (Ow*"/0E) 7 , We find 


4. The velocity maximum. For the function R(é, 7) = 2 s 
R(E, n) = au** + bu*r* cux’ , (4.1) 
where 
a = (Qrkn””)'(n? — 1) {[4n + (2 — 2n + ERIK 
(4.2) 
— (2nk’*)""[8n° + (En — 40°)k? + EL — )k'JE}, 

b = 2n(1 — 1) 7 {T(E, n) — 1 — E(Qrkn””)" {[4y? — (1 — 49? + 2n? + n')KIK 

(4.3) 


1) 
i) 


— (2k’?)"[8n*? — 211 — 2? + 4y° + DR +1 PRE 
ce = —(Qrkn’’’)* {[4 — (2 + nk ]K — (2k) "[8 — (8 + 2n)k? + (1+ a)kJE}. (4.4) 


Near 7T' = 0, the equation of the profile 7’ = constant is 


(p> — 1)(2 — p’*)"” — (3/8)(p° — 1)(2 — p’)'”"T + O(T”). (4.5) 


z= 


The coefficient of 7 is the reciprocal of (2.21), and « = (p’ — 1)(2 — p’)”’”’ is the equa- 


tion of the Blasius-Fuhrmann half-body. 
Expressing y¥* in a similar form for 7 = \* near zero and computing the velocity 


components, we obtain 
wt = (p/2)(8 — 3p’)'” + (3/8)(2 — p*)*(3p° — 2)(8 — 3p) 'T + O(T”). (4.6) 
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From (4.5) and (4.6), 


tm = 6”? — (13/48)6'?T + O(T”), (4.7) 


Winax = 2/3' + O(T”). (4.8) 
The x-intercept x»(7) of the profile is given by 
2(T) = (T — 1)(2 — TT’). (4.9) 


Combining (50) with Fig. 3, we can plot z,, measured from the nose of the half-body as 
a function of 7. 
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SOME IMPLICATIONS OF WORK HARDENING AND IDEAL PLASTICITY* 


BY 
D. C. DRUCKER 


Brown University 


1. Summary. The purpose of this note is to point out the severe restriction imposed on 
possible stress-strain relations by a mathematical formulation of the concepts of work 
hardening and ideal plasticity. Using this condition, an algebraic derivation and ex- 
tension is given of Prager’s extension’ of the Mises plastic potential function. A brief 
discussion is also given of the meaning of stability of plastic deformation as contrasted 
with stability of non-conservative systems in general. 

2. Introduction. A mathematical theory of plasticity must of necessity be based upon 
simplifying assumptions. For example, a good first approximation may be obtained 
by ignoring time effects in most structural metals at room temperature under isothermal 
conditions. A stress-strain diagram, therefore, represents a succession of states of static 
equilibrium. Furthermore, it is an experimental fact that metals generally work harden; 
then each equilibrium state is said to be stable. The phenomenon of the upper yield 
point, exhibited by a few metals, notably mild steel, is an important exception. Another 
exception is the region beyond the ultimate strength on a nominal tensile stress-strain 
plot when necking occurs, Fig. 1. However, this real instability is in terms of force; the 
actual stress continues to increase until failure occurs. 


Oj 
NSTABLE 


STABLE % 


™ Cm 


IDEAL 
PLASTICITY 








€ 


Fic. 1. Plastic deformation. 


3. Work hardening. The meaning of work hardening in simple tension is just that 
stress is a monotonically increasing function of strain, Fig. 1. The plastic deformation 
is then said to be stable. For more general states of stress and paths of loading no such 
simple picture can be drawn. The concept of work hardening, or stability in a restricted 
sense, can be expressed in terms of the work done by an external agency which slowly 
applies an additional set of stresses and then slowly removes them. The original con- 


*Received Feb. 21, 1949. The conclusions presented here were obtained in the course of research 
conducted under Contract N70nr-358 sponsored jointly by the Office of Naval Research and the Bureau 
of Ships. 

1W. Prager, Recent developments in the mathematical theory of plasticity, presented to the 7th Inter- 
nat. Congress for Appl. Mech., London, 1948, J. Appl. Phys. 20, 235-241 (1949). 
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figuration, or state of strain, may or may not be restored. This external agency is to be 
understood as entirely separate and distinct from the agency which causes the existing 
state of stress and which has produced the existing state of strain. 

Work hardening implies that for all such added sets of stresses the material will 
remain in equilibrium and further that 

a) positive work is done by the external agency during the application of the set 


of stresses, and 

b) the net work performed by it over the cycle of application and removal is zero 
or positive. 

[t should be emphasized that the work referred to is not the total work done by all 
of the forees acting, it is only the work done by the added set on the displacements 
which result. Rephrased, work hardening means that useful net energy over and above 
the elastic energy can not be extracted from the material and the system of forces acting 
upon it. 

Consider a unit volume of material in which there is a homogeneous state of stress 
o;; and strain e;; . Suppose an external agency to apply small surface tractions which 
alter the stress at each point by do,; , and produce small strain increments de;; . Next 
suppose the external agency to remove these added tractions, thus releasing the elastic 


strain increments de‘; It then follows from condition (a) that (repeated subscripts 


denote summation) do;;de,;; is positive and from (b) that do; ;(de,; de;;) is either 
positive or zero. Since the plastic strain increments de?; = de;; — de‘; , 
de, (de‘; {. de’ ;) > i. (1) 
do, de; > 0. (2) 


An equivalent manner’ of specifying work hardening which leads to (1) and (2) is 
that as any set of strain increments is increased in ratio from any initial state of stress 
and strain, more and more total and also plastic work per unit volume is required per 
given change in any one of the components of strain. 

4. Stress-strain relations (work hardening). In this section, two basic assumptions 
will be made for a work hardening material. 

1. A loading function exists. At each stage of a plastic deformation a function f(¢;;) 
exists such that further plastic deformation takes place only for f(¢;;) > k. Both f and 
k may depend upon the state of plastic strain and the plastic strain history. 

2. The relation between the infinitesimals of stress and strain is linear. The sum of 
the increments in plastic deformation, obtained separately from each of two sets of 
infinitesimal stress increments do{, and dof; , neither of which constitutes unloading, is 


the same as the plastic deformation resulting from 
do, = dof, + daft ) (of 0o;.)de{, > 0, (0 f/d0,,)dof{ > 0. 


7 . . a . 1 - . 
Following the geometric proof given by Prager, it is necessary to show that any 
do,; representing loading, (0f/d0,:) do,, > 0, may be decomposed into do;, and doi; , 


2D. C. Drucker, The significance of the criterion for additional plastic deformation of metals, presented 
at the 1948 Annual Meeting of the Society of Rheology, J. Colloid Sci. 4, 299-311 (1949). 
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where the arbitrary choice is made that do{, produces no plastic deformation, and 
doi; is proportional to the gradient of f: 


ra] 

OS dat, = 0; 3) 
OO7pR1 

af ) a 

= deli = . f dont 3 deii=r - f , r> 0. (4) 
O07“: OC>%: OC>E: 


This can always be done as (0f/d0,,) dott = r(df/d0,.)(Af/de..) > 0 and can be made 
equal to (0f/Ac,,) do,, by a proper choice of r, 
_ __(0f/Oori)donr_ 
(Of/Oomn)(Of/OFmn) 
The linearity assumption states further that if all do,,; are increased in ratio, all 
de; are increased in the same proportion, so that 


de® ; ~~ Tr 
or 
, af 
de, = gi; —— dow , (5) 
OC7%) 


where the g;; may depend upon stress, strain, and their history but do not involve do,, . 
Using the work hardening condition, we can write Eq. (2) as 
0 < do,,dé; = dot dé; + dotide?; , 
for any do’; such that (df/d0;;) do!; = 0, and for dot} producing de’, . Since the rela- 
tionship do; Cda’, + dot’ also constitutes loading for all values of C, whether positive 
or negative, and produces the same de’; , the work hardening condition becomes 
C dojjdei; + dojjde;; > 0. (6) 


The term Cde’, dé; is the work done on the plastic strain increments by a set of stress 
increments which produce no plastic strain. Its magnitude must be zero, for if it were 
not, a C could be chosen of proper sign and magnitude to violate (6). 


Therefore, 


det, = debits Zt iy wt 
or 

do! gi; = 9, (7) 
where 

0 

gt do!; = 0. (3) 


Conditions (7) must be satisfied for all do‘; satisfying (3) so that 


_ of 
gis =G i (8) 
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where G is a scalar which may depend upon stress, strain, and history. The only per- 
missible stress-strain relation consistent with the two basic assumptions made at the 
beginning of this section is, therefore, 


‘ 9 f 
de’; G = f z } ds, : (9) 
O60 ;; 00%; 


It is of interest to note that the work hardening condition (a) automatically insures 
that unique increments of stress in the interior of a body result from given increments 
in surface traction’, or, more generally, that the history of the state of stress at each 
interior point is uniquely determined by the history of the surface loading. 

5. Ideal plasticity. Ideal plasticity in simple tension is indicated by a horizontal line 
on the stress-strain diagram, Fig. 1. The stress remains constant as the strain increases. 
For more general paths of loading, ideal plasticity means that the work done by an 
external agency which slowly applies and then removes a set of stresses is zero over an 
equilibrium cycle. Equation (2), for a work hardening material, is replaced by 

do, ,de{; = 0 (10) 
for an ideally plastic one. The equality in Eq. (2) applies only when de?; = 0 so that 
the difference is by no means trivial. 

6. Stress-strain relations (ideal plasticity). The further assumption of ideal plasticity 
is that f(o,;), a function of stress alone, exists such that plastic deformation takes place 
without limit when f(c,;) = k and that the material is elastic for f(o,;) < k. During 
flow, therefore, 

Of 


005; 


de.; = 0. (1 1) 


7 


As all do;; satisfying Eq. (11) must also satisfy Kq. (10), 


0 
, of (12) 


Pp _— 
dé; = . 


005; 
where \ is a scalar. The form of f may be as anisotropic as desired and of any degree in 


the components of stress. 
7. Stability of non-conservative systems. When a metal work hardens, the plastic 





SA 
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Fic. 2. A simplified model of plastic flow. The 

coefficient of friction is variable for a work 

hardening material. Here do;; de?;; > O is 
required. 


deformation is termed stable. The question arises as to whether this is stability in the 
conventional sense. Figures 2 and 3 show that work hardening conditions (a) and (b) do 


imply stability, while Figs. 4 and 5 show that the converse is not true. 
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Suppose, in Fig. 2, that P is considerably less than the maximum permissible fric- 
tional force F (conservative system case). Application of a small positive or negative 
increment of force AP will produce a displacement of the end of the spring in the direc- 
tion of AP. The external agency adding AP will do positive work in each case and will 





. | WAS 


rrr eS pie it iin 44444 474?, 




















F1a. 3. A model in which the elastic strain is not 
recoverable directly. Here do;; (det;; + de®;;) > 
0 is required. 


do zero net work over a cycle of application and removal of AP. For a conservative 
system, work hardening condition (a) is equivalent to a minimum energy principle. 
Although not useful for such a system, condition (b) is always satisfied, as the net work 
is Zero, 

Now let P equal the maximum permissible value of F for the position of the block, 
and further suppose that the coefficient of friction will increase as the block moves to 
the left. Under a decrement of force, the system acts in a conservative manner. How- 
ever, an increment AP will stretch the spring and also move the block; for both reasons, 
AP will do positive work. Remova! of AP releases the elastic energy increment, but 
the block remains in place. Consequently, the net work done by the external agency 
over the cycle is positive. The system is obviously stable under disturbances AP; an 
input of energy is required to displace the block. If the coefficient of friction decreases 
as the block moves to the left, negative AP is required to maintain equilibrium when 
the motion occurs. The work done by AP is negative, and the system is clearly unstable. 
In this case, an equilibrium cycle is not possible unless the friction coefficient should 
increase again, as the motion proceeds. 

The system of Fig. 3 gives comparable but somewhat different results. Suppose P 
to be equal to the sum of the maximum frictional force F for the position of the block 
and the spring force S. If F remains fixed or increases with displacement to the left, 
an increment AP in P does positive work. On removal of AP, the elastic energy increment 
remains in the spring. Work hardening conditions (a) and (b) are satisfied, and the 
system is clearly stable under possible disturbances AP; an input of energy is required 
for motion. On the other hand, if F should decrease more than S increases, the system 
is unstable. Work done by the external agency during the application of the additional 
force required for equilibrium and over the complete cycle is negative. 

Figures 2 and 3 illustrate the sufficiency of conditions (a) and (b) for stability of 
the usual type; no displacement from the equilibrium position without an input of 
work. They are: one-dimensional examples, however, and the analogy between P and 
stress and between displacement and strain is really not so close. 

Figure 4 brings out the dissimilarity and also the more general problem of stability. 
If V remains constant, Fig. 4 reduces to Fig. 2. However, if P = uN and V is increased 
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by AV, P must be decreased by nAV to maintain equilibrium. Positive work is done 
in the application of AV and —yAV. If these added forces are removed in ratio, the 
total work over the cycle is zero. If AV is removed first, the result is the same; but if 
—yAV, the decrement in P, is taken away first, the block will displace to the left and 


Vv 





BAW 
YW 
——— F< J N 
N=W-V 














Fic. 4. Work can be extracted. 


negative work will be done by the external agency over the complete cycle of application 
and removal of the additional forces. Repeated cycles will cause continued displacement 
of the block if the friction coefficient remains constant, but will produce no additional 
motion, after the first cycle, if 1 increases. From one point of view the system is un- 


p s 


A i "SA 


(a) (b) 


Fic. 5. Hydrostatic pressure and shear. 





























stable under the disturbances because the block will move without any external agency 
doing work on it. More realistically and usefully, if the block moves only infinitesimally 
for repeated infinitesimal disturbances (u increasing), the system would be considered 
stable. The case of constant ~ would be termed neutral, and only if finite work could be 





STABLE 








Fic. 6. Stability of non-conservative systems. 


extracted for infinitesimal disturbance (u decreasing) would the term unstable be justified 
in general for a non-conservative system. 
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Another simple example of these distinctions is furnished by a block on a surface 
inclined at the friction angle, Fig. 6. Assuming a constant coefficient of friction, if the 
surface is a plane, infinitesimal disturbing forces applied to the block normal to the 
plane will cause the block to move and continue to move down the plane for any finite 
distance. During this downward motion the block cannot do a finite amount of work 
on an external system, as the work done by the gravitational force is dissipated by the 
frictional work. The term neutral seems appropriate. If the surface is curved so that 
the slope becomes steeper as the block moves, the block will accelerate. Work can, 
therefore, be done on an external system; for example, a spring may be compressed. 
The system is unstable. If, on the other hand, the curvature of the surface diminishes 
the slope, then an infinitesimal disturbance will produce an infinitesimal motion. Repe- 
tition of the disturbance will cause no additional motion; the system is clearly stable. 

The stability of plastic deformation implied by work hardening is stability in the 
strictest sense. Even the infinitesimal motion for Fig. 4 with uw increasing is not per- 
missible. Figure 5 emphasizes this feature for plastic deformation and is a truly critical 
illustration of the concepts advanced. Suppose Fig. 5a to represent a portion of a work 
hardening body subjected to a homogeneous state of stress composed of a hydrostatic 
pressure p and a shearing stress r. Consider the external agency to apply and then re- 
move a hydrostatic tension s and shearing decrement ¢. As Bridgman has shown,’ a 
high hydrostatic pressure increases the flow stress of metals. For simplicity, assume the 
yield stress in simple shear at some stage of plastic deformation to be given by 


Ty, = To + kp, (13) 


where 7, is some initial value and k is a positive constant. At p — s the yield value will 
be to + kp — ks. Now let + = 7, and t = ks. The states of stress given by Figs. 5a and 
also 5a plus 5b are then exactly on the yield or loading surface. If t is removed while s 
is maintained, the body will deform plastically. Let the added hydrostatic tension be 
removed next. If the assumption is made that the flow pattern is unaffected by the 
hydrostatic pressure and is, therefore, simple shear strain, work will be extracted from 
the system in the cycle. The same result can be obtained directly by taking s,; = s 
but ¢; < ks. When s, and ¢, are applied simultaneously instead of s and t, negative 
plastic work is done by ¢, , zero plastie work by s, ; do;; dé; < 0. 

This contradiction to Eq. (2) shows that such a material is analogous to the friction 
example, Fig. 4, and is not truly work hardening despite the fact that tensile or shear 
stress-strain curves would look entirely proper. However, the conclusion is not that work 
hardening metals do not exist but is instead that real metals which do exhibit dependence 
on the mean normal stress also increase in volume as the plastic deformation proceeds* 
and that the plastic work done by s or s, is more than that done on ¢ or ¢, in most, and 
possibly all, cases of engineering interest. 

Taking the loading function f in Eq. (9) as 


f=4+d’ +aJ, (14) 





P. W. Bridgman, Effects of high shearing stress combined with high hydrostatic pressure, Phys. Rev. 
48, 825-847 (1935). 

‘A. H. Stang, M. Greenspan, and 8S. B. Newman, Poisson's ratio of some structural alloys for large 
strains, J. Research Nat. Bur. Standards 37, 211-221 (1946). 
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where a is a constant, J, 8;;8;;/2, J: = Orn » 83; o;; — J,6;,;/3, illustrates but does 
not, of course, prove this point. Loading, or further plastic deformation, results if 

of 1/2 
: do; ; Aids )+a dJ, ~ 0, (15) 


O05; 


which always insures positive 


do 


. 7 of of 
de’, G do; ; — (. , dou). 
06>; OO>%: 


ij 
The plastic stress-strain relation is 

p Y 85; ¢ . 
dei; = G\ s7ia + 26;;) af, (16) 


and the plastic dilation is, therefore, 
det, = 3aG df. (17) 
For the loading shown in Fig. 5, the plastic increment in the engineering shearing 
strain is 
dy’ = 2de, = G df, (18) 
so that 
dey, /dy’ = 3a, (19) 


where 3a is the small slope of the straight line relating the yield stress in shear to the 
mean normal pressure, k in Kq. (13). 

The preceding example shows that consistency is obtained and uniqueness of solution 
assured if the small plastic volume increase which does occur is taken into account. It 
does not, however, prove that some material system which can do work on an external 
agent in a closed force cycle, does not, or can not, exist. The argument is only that such 
a material is not work hardening. Moreover, it is found experimentally in the plastic 
deformation of most metals in the range of engineering interest that the only energy 
available is in the form of elastic energy, and that no additional amount of work can 
be obtained from the system of forces acting on a plastically deformed metal which is 
in a state of stable equilibrium. No work can be extracted in a force cycle, since the 
elastic energy is the same at the beginning and the end of the cycle. 

8. Conclusion. The usual basic assumptions of the theory of plasticity when coupled 
with a mathematical formulation of the concepts of work hardening and ideal plasticity 
require stress-strain relations to be of the Mises-Prager form, de?; = A(df/d0;;). For 
a work hardening material, \ = G(0f/d0,,) do,, . G and f may depend upon stress, 


plastic strain, and the history of plastic strain. 
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—NOTES— 


ON THE CONVERGENCE OF MATRIX ITERATION PROCESSES* 
By ROBERT PLUNKETT (The Rice Institute) 


A recent paper [l1]t has presented a formal proof of an iteration method for the 
solution of non-homogeneous (Fredholm) integral equations which is directly applicable 
to non-homogeneous matrix equations. Various special cases of this method have been 
used for matrix equations (ef. [2], [3]) but no general statement has yet appeared. The 
proof in the case of matrices is much simpler and the conditions for convergence are 
less stringent than in the case of integral equations. 

The equation under consideration is 


U =U + Qu, (1) 


where u is an unknown column vector of n elements, up is a known column vector and 
Q is a square matrix of n* elements. The formal solution is obviously 


u= (I — Q)'w, (la) 


but the reciprocal matrix is not always easily determined if n is very large. Associated 
with Eq. (1) is the homogeneous equation 


A,u, = Qu,, (2) 


in which A, are the latent roots and u, , the corresponding latent vectors. 
The usual iteration method is to let 


UV, = Uo + Qua-s (3) 


where Up u + w is the first approximation to the answer and », is the Ath iterated 
approximation. For rapid convergence, an attempt is made to make w as small as 
possible. Then 


Vv, = Uy + Qu = u + Qu, 
(3a) 
v», =ut Q’w. 


For a matrix Q with n distinct, non-repeated roots, Eq. (3a) is easily expanded into 
a simpler form. When all the roots are distinct, all the latent vectors are linearly inde- 
pendent, and w may be expressed uniquely as a linear combination of them: 
n 


w= bs As 


r=1 


where A, are scalars, and it is obvious that 


Q’w = > r*A,u, , (3b) 


r=l1 


which must converge to zero for all | A, | < 1. 


*Received Feb. 21, 1949. 
tNumbers in brackets refer to the bibliography at the end of this paper. 
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If the roots are not all distinct, let us assume that there is one \, which is repeated s 
times, and that | JA, — Q| is q times degenerate, i.e., that there are q linearly inde- 
pendent vectors associated with this root. If g = s, the equation may be treated as in 
Eq. (3b). If g # s, the confluent form of Sylvester’s Theorem is applicable and states 
(ef. [4]) that 

P(Q) = >> T(a.), (3c) 
the summation to be taken over all distinct values of \, , where P(Q) is any polynomial 
in Q and 


P'(d.)Z.~s(r,) 
1! 


l' (0 
Si) ws = E “ )] ; 
UL ay* \GAA) ni, 


Now P*'’(\,) is the 7th derivative of P with respect to \ at A, , F(A) is the adjoint of 
(IX — Q), and 


pis a ( J e 
4 I A, )Zo(A,) 


T(A,) P(A,)Z,-1(\,) + be sine df 


’ 


(s — 1)! 


A <-D —%,,)0 —'ad)* A~&) 

i.e., the products of the differences involving all the roots except the s repeated roots 

under consideration, but including the repeated values of all other repeated roots. Letting 
P(Q) = QY = > 7,(.), 


then 


T,(X.) = AZe-1(A.) + tw # ae ae: 


(h - Serre | —~s ¢ h—#+2 
Ae vce } + 1)/ h oa 3X Zr ’ 


(s — 1)! 
and for large h it is apparent that 
Q**'w = bh Ty41(A,)0 = } A,T,(A,) ¥, 
since the matrices Z do not depend on h. The validity of this last expression is increased 
for multiple degeneracy of (JA, — Q); if the degeneracy is g, F(A,) and all of its de- 
rivatives up to the (¢ — 2)th are zero, which means that 


Z(A.) = 0 (« < q — 2) 


or that the last g — 1 terms in T(A,) vanish. This is to be expected, since for degeneracy 
s, Eq. (3b) holds exactly. Then for large h, Eq. (3a) may be written 
Q’w = Be d’B,w, (3d) 


where B, is a matrix involving the Z;(A,) and is independent of h in the limit, and the 
summation is taken over all distinct values of \, . Thus it may be seen that Eq. (3a) 
converges to u for large h in all cases, if all | A, | < 1. 
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The proposed solution is to let 


M% = O%-1 + (1 — 6) [uo + Qv,-1], (4) 
where 
Vv, = Wo + [uo + Quvoj(l — 4) 
=ut+ [al + (1 — AQ)w 
or 


v», = ut [al + (1 — 0)Q]*w. (4a) 


Comparing this with Eq. (3a), it may be seen that the convergence depends on the roots 
of [a7 + (1 — @)Q] which are [6 + (1 — 86)A,]. Thus if all 


|ea+(1— @)A,| <1, (4b) 


Eq. (4a) will converge to wu. 

There are cases in which it is impossible to find a number @ which will satisfy Eq. 
(4b) for all A, , without manipulating the Eqs. (1) so as to change Q and thus the roots. 
For example, if there are two of the roots which are real, such that 0 < A, < 1 <\X,, 
then for A, , @ < 1 and for \, , @ > 1 which is impossible. But if all A, are real and less 
than 1, where A, is the largest negative root, (| A, | — 1)/(/A,| +1) < 6 < 1. Similarly, 
for all A, real and greater than 1, 1 < @ < (A, + 1)/(A, — 1), where A, is the largest 
positive root. While the exact value of @ within the limits is not important, it may be 
seen that all values of | @ + (1 — @)A,| should be as small as possible for best con- 
vergence. For example, in one actual case [2] there was one root of —2, two near —1 
and several negative and very small, @ was taken to be 1/2, the largest value of 

@ + (1 — @)A, | was about 1/2 and five iterations reduced the error to (1/2)° or about 
3%. 

If all the \, are real, positive and less than one, this method is not necessary for con- 
vergence, but hastens it if any A, are close to one. 

Equation (4) is equivalent to a Toeplitz summation on the vectors v, as specified 
by Eq. (3), with the summation matrix (ef. [5]) 


tia _— 6)" 6’, j < i, 


= Q, j>t. 


a;; 
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A GEOMETRICAL INTERPRETATION OF THE RELAXATION METHOD* 
By O. BOTTEMA (Technical University, Delft) 


1, 2, --+ , m) be given constants and consider the equations 


9 &) 


>, a;;7; — B; = 0, (i = 1,2, --- ,n). (1) 


As n increases, the solution by means of determinants becomes burdensome. Then 
the relaxation method provides a set of easy steps by which the solution is approached. 

Synge** has given a geometrical description which illustrates the relaxation process. 
He confines himself, however, to the special case where the coefficients are such that 

a;; and where, moreover, the quadratic form bs @;;X;2; 18 positive definite. In 
interpretation a set of homothetic ellipsoids in n-space plays a fundamental part. 
The purpose of the present note is to give a simple geometrical interpretation for 


a 


7 


his 


the general case in which the matrix a;,; is not necessarily symmetric. 
If we substitute in the left member of (1) the arbitrary set of values y; , we get 


> ay; — B; = RY, (2) 
i=1 


where FR; are the so-called residues. It is the task of the relaxation process to construct 
a sequence of successively “better” sets of values so that the residues tend to zero. 


Following Southwell, we do this by correcting one of the quantities y, say y, , while 
the others are left unaltered. If y, is replaced by y, + d, , we get the new residues 
Ri? = RY? + aud, . (3) 


Now let us regard R; as the rectangular Cartesian coordinates of a point P in Euclidean 
n-space. The purpose of the process is to approach the origin O of the coordinate system. 
$y means of (3) the arbitrarily chosen point P, is replaced by P, . It is obvious that the 
direction of the line 1 = P,P, is determined by the given coefficients; this direction is 
specified by a, (¢ = 1, 2, «++ , ”). 

The procedure which we have to follow is not completely defined, but it seems 
reasonable to move P from P, in the given direction in such a way that its distance to 
the origin becomes as small as possible. In this case d, has to be chosen so that the 
quadratic function 


dD ahs + 2d, D aul” + D (Ri) 


is 2 minimum. Hence d,; = (—>oan.Rs”)/(d 0%), and this is indeed the method fre- 
quently used; it is in accordance with the principal of least squares. But it is clear that 
P, is the orthogonal projection of O on the line 1; in other words, P, 7s the projection of 
P, on the (n — 1)-dimensional space through O which is perpendicular to the line l. Since 
y, has now been corrected, we take another of the unknown quantities and proceed in 


the same way. 
It is thus obvious what kind of construction in our geometrical representation is 


*Received May 27, 1949. 
**J_ L. Synge, A geometrical interpretation of the relaxation method, Q. Appl. Math. 2, 87-89 (1944). 
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the analogue of the procedure. We consider the (n — 1)-dimensional spaces V,,; 
(i = 1, 2, --- , n) through O, whose equations are bs a,;,2, = 0. When we take arbitrary 
values for the unknown quantities and correct them in the order a,, , %, , **: , the 
construction runs as follows. An arbitrary point P, is chosen in n-space; P, is the or- 
thogonal projection of P» on V,, ; P2 is the projection of P, on V,, , etc. And so the pro- 
cedure is illustrated in a simple geometrical way. We add some supplementary remarks. 
The point P,, (m > 0) lies always in one of the spaces V, . These spaces are linearly 
independent if | a;, | # 0. At each step (after the first) one of the V,, is projected by 
parallel projection onto the following one. In this way, an affine correspondence is 
established between the two successive V,, , the modulus of the affinity being cos a, 
where @ is the angle between them. Thus the convergence of the procedure can easily 
be proved, provided that the corrections take place in a fixed cyclic order 2, , 2%, ,°°* , 


x,, Where k, , kp, +++ , k, is a permutation of 1, 2, --- , n. 
If two successive V,, are perpendicular to each other, the projection of the first 
onto the second coincides with their (n — 2)-dimensional space of intersection. It 


follows, therefore, that if all the V,, are mutually perpendicular (that is, if the matrix 
|| ay, || is orthogonal), the point P, lies on V,, , P2 on the intersection of V,, and V,, , 
P, on (V,, , Ve, , Ve,) and so on; hence P, coincides with O. In this case the procedure 


ends automatically after n steps. 


A SIMPLIFIED METHOD OF DIFFERENTIATING AND EVALUATING 
FUNCTIONS REPRESENTED BY FOURIER SERIES* 


By A. M. WINSLOW (University of Washington) 


1. Introduction. This paper shows how to eliminate the difficulties caused by discon- 
tinuities of Fourier sine series at the ends of the interval of periodicity. 

Applications of Fourier series to exact solutions of problems in mathematical physics 
involve the following essential considerations. In an interval —a S z S a, it is assumed 
that a function f(x) and its successive derivatives up to some finite order f‘” (x) all comply 
with sufficient conditions of continuity, bounded variation, differentiability and in- 
tegrability. They thus permit representation by Fourier series, which can be differentiated 
to give the derivative of next higher order, and integrated to give an expression for the 
derivative of next lower order. 

When f(x) is an z-odd function, and f(a) # f(—a) + 0, particular difficulties are 
encountered. The corresponding Fourier sine series is discontinuous at x = +a and 
does not conveniently represent the values of f(a) and f(—a). In addition, the derivative 
f'(x) is represented by a complicated Fourier series which is not readily evaluated at 


x = a. Thus 
f(z) = do bi sin B,2, (1) 
1 


in which 8, = n/a. The expression for the derivative is 


*Received March 28, 1949. 
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fa) =F (f@ - (-a) + ¥ a [f(a) — f(—a)] + uss, cos Br. (2) 


2. Modification of Fourier sine series. These difficulties of evaluation and differentia- 
tion are alleviated by rewriting Eq. (1) in the form 


f(q) =cxr + (> b/ sin B,7 — cr), (3) 


1 


where c = f(a)/a. The last term of Eq. (3) may be written 
t= z. b/’ sin 6,2; 
1 


then 


Zz b/ sin B,t% — cz = > b, sin 6,2, (4) 
1 1 


in which b, = bi — bj’. 

In Eq. (4), be, b, sin 6,2 represents a function having the value zero at x = ka. 
It consequently permits termwise differentiation. 

Thus the original equation (1) is now 


f(x) = cx + ps b, sin B,2. (la) 
1 
At x = a, f(x) = -:ca, respectively. Also the derivative f’(x), correctly obtained by 
termwise differentiation, is 
f(x) =ct+ z a, COs B,2, (2a) 
1 


where a, = D,8, . 

In this manner, with complete generality, Eqs. (1) and (2) can always be rewritten 
in the form of Eqs. (la) and (2a). 

Furthermore, in the case of integration of an 2-even function f’(x) expressed by 
Eq. (2a), termwise integration is a legitimate process. If it is known that the indefinite 
integral f(x) is an z-odd function, there is no additional constant of integration, and 
f(z) is expressed by Eq. (la), in which b, = a,/8, . Also at x = sta, f(x) = ca, re- 
spectively. This follows from the fact that termwise differentiation of Eq. (la) gives 
Eq. (2a). Therefore, >>> b, sin 8,7, complying with the necessary condition of termwise 
differentiation, must represent a function having the value zero at x =-La. 

3. Derivatives of any finite order. More general results for a mixed function are 
obtained by successive integrations of f‘”(x) with introduction of constants of inte- 


gration. Thus 
f(a) =c’+ >. (af cos B,x + bf sin B,x) (5) 
1 
after m successive integrations gives an expression of the form 


f(z) = DSe2" + Dd (a, cos Bx + b, sin 8,2). (6) 
0 1 
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Conversely, derivatives up to the order f‘”’ (x) may be obtained from Eq. (6) by term- 
wise differentiation. In the expressions for f(z) and any derivative up to f‘"~’’(x), any 
sine series represents a function having the value zero at x = -ta, thus facilitating 
evaluation at these points. In some solutions Eq. (5) is a convenient expression for 
f° (x). If desired, Eq. (5) can be modified by adding a term c’’x combined with >>? 
b/ sin 8,2 so that the sine series in this equation also represents a function having the 


value zero at x = +a. Then Eq. (6), modified by this additional term, is 
m+1 @ 
f(x) = Zz c,2 + 2. (a, cos 6,2 + b, sin 8,2). (6a) 
0 1 


Equations (6) and (6a) are special cases of Borel’s theorem’ which apply to derivatives 
of any finite order up to f*”’(z). 

The preceding methods of derivation, stated in terms of a single variable zx, also 
apply completely to a three-dimensional member bounded by parallel planes x = +a 
with assigned boundary conditions. In this case f(x) is replaced by f(z, y, z) and the 
derivatives are written as partial derivatives with respect to xz, f°” (x) becoming d”f/dz”. 
Also in Eqs. (la), (2a), (6), (6a), every coefficient a, , b, , c, , c, While independent of z, 
is a function of y and z. For example, in Eq. (la), 6, is b,(y, z), and ¢ is c(y, z) = 
f(a, y, 2)/a. Thus f(x, y, z) and each of its x-derivatives of finite order can be expressed 
as the sum of a finite power series and a Fourier series without discontinuity at « = +a; 
the expressions can, therefore, be evaluated definitely at these boundaries. 


1. Borel, Lecons sur les fonctions de variables réelles, Gauthier-Villars, Paris, 1905, p. 68. 


A TAPERED LINE TERMINATION AT MICROWAVES* 
By GEORGE J. CLEMENS (City College of New York) 


1. Introduction. In the field of ultra-high frequency, one method of power trans- 
mission is by means of concentric metallic conductors called a co-axial transmission 
line. The maximum amount of power can be delivered over these lines if the impedance 
of the load is equal to the characteristic impedance’ of the line. In addition, in various 
measurements at ultra-high frequency it is essential to have a matched termination 
over a broad band of frequencies. This problem was approached experimentally and 
led to the tapered line termination as shown in Fig. 1. 

To the left of z = 0, the co-axial line has an inner metallic conductor of radius a. From 


x = 0 tox = L, the inner conductor is a glass tube coated with a thin metallic film of 
resistive material.” To the left of z = 0, the outer conductor has a radius b and from 
x = 0 tox = L, the outer conductor has a linear taper down to a radius c. Also at 
x = L, the inner and outer conductors are short circuited. 


*Received Dec. 2, 1948. The results and general method of this paper conform to a section of the 
author’s doctoral thesis, the research for which was carried on under the direction of Dr. Ernst Weber at 
the Polytechnic Institute of Brooklyn. 

1J. C. Slater, Microwave transmission, McGraw Hill, 1942, pp. 71-74. 

2Publication P. B. 6588 (U. S. Commerce Dept.), 1945. 
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2. Differential equation of the problem. The voltage and current distribution in a 

co-axial line can be written as a pair of differential equations® 
dV/dx = IZ, dI/dx = —VY, (1) 


where V and 7 are the voltage and current at any point x on the transmission line, 
and Z and Y are the series impedance and shunt admittance per unit length at any 






































2a 
Fie. 1. 
point x. Differentiating each of Eqs. (1) and combining, we obtain a pair of second 
order differential equations 


d?I/dx? — (Y’/Y) dI/dx - ZYI =0, 
(2) 


I 
S 


d’V/dx? — (Z'/Z) dV/dx — ZYV 
where Y’ and Z’ are the derivatives of Y and Z with respect to x. In general, with air 
filling the space between the conductors, we find that 

Z = R + wwk, log (r,/r,), 
| (3) 
. Ww v 
Y = 9g + i——- , 
: log (To r;) 
where R = resistance per unit length, g = conductance per unit length, w = 2zf, f = 
frequency, k, = 2re, = 10°°/18 farads/meter, k, = uo/2r = 2 X 10°’ henries/meter, 
ro = radius of the outer conductor between x = 0 and x = L, r; = radius of the inner 


conductor between x = 0 and x = L. 
The actual linear tapered line shown in Fig. 1 is now approximated by an exponential 


line in which the variation of the outer conductor is given by 
Tf) = a(b/a) exp (—kz). (4) 
We then find that (assuming g = 0, R constant) 
Z=R+ Z, exp (—kz), Zo 


iwk, log (b/a), 
(5) 
Y = Y, exp (kz), Y, = twk,/log (b/a), 
where Z, and Y, are, respectively, the series impedance and shunt admittance per unit 
length of a lossless uniform line (R = g = OQ) having a constant outer radius b and 
constant inner radius a. The differential equation (2) for J will then become 


@I/dx? — k dI/dx — [YZ + RY, exp (kx)|I = 0. (6) 





3J. C. Slater, loc. cit., pp. 60-71. 
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If we let 
s = [—4RY, exp (kx)/k’]'” (7) 
and set J(x) = sW(s), we can transform Eq. (2) into 
d’W/ds* + (1/s)dW/ds + [1 — (4¥0Zo/k*)(1/s")|W = 0. (8) 
This differential equation for W(s) is the Bessel equation of order 
vio = £(4Z,Yo/k’)'”. (9) 
Hence, we may write 
I(x) = exp (ka/2)[MJ,,(s) + NJ,,(s)]. (10) 


It should be pointed out that the substitution, Eq. (7), which resulted in the differential 
equation (8) was obtained only after a series solution of an integral equation was found. 
3. Integral equation. The method of solution for the lossy exponential line by means 
of an integral equation was suggested by a “perturbation method” of Schelkunoff.* 
If we consider a lossless (R = g = 0) exponential line, then 


Y = Y, exp (kz), Z = Z, exp (—kza), (11) 
and Eq. (2) becomes 
d’I/dx® — k dI/dx — ZY I = 0. (12) 
Since k and Z)Y, are constants, the solution for I(x) is 


I(x) = A exp (m, x) + B exp (m,2), 


where 
m, = 1/2[k + (kh + 4T%)'”’], 
ae ead Yo MQ - (13) 
m, = 1/2[k — (kh? + 473)"”’), 
Using Eq. (1), we find that 
V(x) = —Yo' exp (—ka)[m,A exp (m,x) + mB exp (m;2)]. 


We now assume that at x = «, there is a point discontinuity of R de ohms in series 


with the line. To the left of z = e 


I,(x) = A exp (m,x) + B exp (m2), 


(14) 
V(x) = —Y,' exp (—ka)[Am, exp (mx) + Bm, exp (m2)]. 
To the right of xz = e 
I,(xz) = C exp (m,x) + D exp (m2), 
(15) 


V(x) = —Y,' exp (—ka)[m,C exp (m,x) + m2D exp (m2)). 








‘S. A. Schelkunoff, Electromagnetic waves, D. Van Nostrand Co., 1948, ch. 4. 
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At x = «, the boundary conditions for V and J become 


Ine — 0) = Ife + 0), 
(16) 
Viole — 0) = Vile + 0) + VO), 


where V(e) = J(e)R de. We may now determine the constants C and D in terms of A 
and B. If we assume an infinite distribution of discontinuities along the line, we obtain 
the following integral equation for J: 


I(x) = Ine) — RYo(m, — m)" | I(© exp (Ke) 
“0 (17) 
-{exp [m,(a — &)] — exp [m.(x — 6)]} de. 


In Eq. (17), R& is constant. 
Using Picard’s method of successive substitution” and rearranging terms, we find 
the infinite series 


Te P exp (kx) | _P" exp (2kx) _ 
Iz) = E exp (ma)| 1 ¥ 2m, k + 4m,k*(k + 2m,) 
P* exp (3kz) 
a c.f. ee 8 
> 12m,k°(k + 2m,)(2k + 2m,) +r | (18) 


hae oil P exp (kz) , _P* exp (2ka) _ | 
+ F exp (ma)| + 2mk + Am2k*(k + 2mz,) . 


where E and F are constants and P = RY, . With a further rearrangement of the series 
and a change of variables we obtain the Bessel solution of Eq. (10). 

4. Comparison of exponential and linear tapered lines. In order to compare the expo- 
nential and the linear tapered lines, we must find the input impedance to the exponential 
line. This is obtained by using Eq. (1) and finding the quotient V(x)/I(x) at x = 0. 
From Eq. (10) we find that 

V(x) = (—1/Y){M exp (kx/2)[(k/2)J,,(s) + J?,(s)(ds/dz)] 
(19) 
+ N exp (kx/2)[(k/2)J,,(s) + J/,(s)(ds/dz)]}. 
Applying the recurrence formulas for Bessel functions,” and using the relation (Eq. (7)) 
ds/dx = ks/2, we find that 
V(x) = [—k exp (kx/2)/2Y]{M[(v, + 1)J,,(8) — sJ,,+:(8)] 
(20) 
+ N[{(v. + 1)J,,(s) — sJ,,41(8)]}. 
Since the transmission line was shorted at x = L, we find that V(L) = 0, and from 
Eq. (20) we have that 





5B. G. Keller, Mathematics of modern engineering, vol. 2, John Wiley, 1942. 
6N. W. McLachlan, Bessel functions for engineers, Oxford Univ. Press, 1941, p. 158. 
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M/N = —[(v, + 1)J,,(sr) = SrJ,,+1(81)] ‘T(v, + 1)J,,(sz) _ SrJ,,41(82)], 
(21) 
= [—4RY,k™ exp (kL)]"” 


Then the input impedance becomes 


Zi. = V(0)/I1(0) = T(1 — U)/(1 — W), 





where 
po = gb [act Mal) — eld] 
7 2Y, J,,(80) 4 
y [(v, a. Liz. (s,;) = oF J. wise ICA —- I). - (sp) —~ Sod » +1(8) ] 
T — e -~ cs 2 — 43 A 8 1 1 
(vy, + 1)J,,(s1) — etal t)\M(2 + 1)J,,(80) — 8od,+1(80)] ’ (22) 
WT = ws Ad ue. —_SL J), +1(81)] J», (80) 


i+ ‘Du on SL , = i(sz)| J ,,(8o) 


A measure of the input impedance of a termination is the voltage standing wave 
ratio’ (VSWR) produced on the concentric co-axial line to which the impedance is 
attached. The VSWR is the ratio of the maximum voltage to the minimum voltage 
between the inner and outer conductors. If the impedance given by Eq. (22) is the termi- 
nating impedance on a transmission line having a characteristic impedance Z, , then 
as shown by Slater,*® the VSWR produced by this termination is given by 


VSWR = (1+ |K))/—-|K)), 
(23) 
K = (Z, — Zin)/(Z. + Zin). 


In general, Z, , Z;, , and K are complex quantities. In order to check the approximating 
exponential termination against the linear tapered termination, it would be necessary 
to calculate values for Bessel functions of complex argument and indices. At present 
there are no such tables available. Asymptotic expansions in the literature as well as 
those mentioned in Watson’s treatise on Bessel functions” are either not applicable or 
involve series which are no more rapidly convergent than the actual series themselves. 
Since we wish to determine how the input VSWR varies with the wave length of the 
impressed frequency, it is necessary to modify Eq. (22) in terms of A, the wave length. 


Since 
rs: = 1 = —(2m)* f’uo€o = —(2mr)* f'v” = —(2r/d)’, 


we find from Eqs. (7) and (9) that 


> 


v1.2 = +1[(4r/kd)’ — 1)”, 
(24) 
(2/k) exp (ir/4)[22R exp (kx)/dZ,]'”. 


~” 
II 


7J. C. Slater, loc. cit., pp. 48-54 and King, Mimno, Wing, Transmission lines, antennas, and wave 
guides, McGraw-Hill, 1945, ch. I. 
8G. N. Watson, T'reatise on the theory of Bessel functions, Cambridge Press, 1945. 
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The numerical check is made against the linear tapered line having a length of 5.5 
centimeters and a total resistance of 59 ohms. This termination is connected to a co- 
axial line having a characteristic impedance of 46.4 ohms. (a = 0.406”, b = 0.1875’). 
If & is equal to 1.3/L, a comparison of the radii of the tapered outer conductors is given 
in Table 1. 

TABLE 1 


Rapius oF OuTER CONDUCTOR 


z/L Linear Taper Exponential Taper 
0.0 . 4067 in. . 4067 in. 
0.1 .385 .370 

0.2 . 364 . 342 
0.3 . 342 .317 
0.4 .321 .301 

0.5 . 300 . 282 

0.6 .278 . 269 

0.7 . 257 . 256 

0.8 236 246 

0.9 .214 . 238 

1.0 .193 , 232 


If in Eq. (24) we use the following values: 
k 


R = 59/0.055 = 1072.7 ohms/meter, 


1.3/0.055 = 23.636, 


Z. = 46.4 ohms, 
we will find the following comparison between the measured VSWR for the actual linear 
tapered line and the calculated VSWR for the exponential line. 
TABLE 2 
Measured VSWR Calculated VSWR 


d (meters) f (megacycles) Linear Tapered Line Exponential Line 
0.303 990 ee 1.26 
0.2796 1073 1.12 1.29 
0.2545 1179 1.12 1.31 
0.2389 1256 1.12 1.33 
0.2108 1423 1.14 1.33 
0.1993 1505  Y 1.27 
0.1200 2500 1.13 1.26 
0.1111 2700 1.12 1.19 
0.1001 2997 1.12 1.18 
0.0899 3337 ee 1.16 
0.0805 3728 1.08 1.13 
0.0756 3968 1.09 1.14 


5. Discussion of results and design factors. The results of Table 2 indicate a rela- 
tively flat response for the VSWR over a broad band of frequencies. The differences 
between the calculated and the measured values are due to the following reasons: 

1. Since the exponential line cannot duplicate exactly the linear tapered line, we 
will of necessity have a slight discrepancy. Now it is much more difficult to machine a 
tube with an exponential variation for its inner radius than to machine a straight taper. 
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Moreover, since the exponential taper indicates a conservative result, it would be more 
practical to use the linear taper. 

2. An assumption was made in setting up the original differential equation. In a 
rigorous solution of the electromagnetic fields existing between a pair of concentric 
cylinders,’ only the principal mode will be transmitted if the wave length of the im- 
pressed field is greater than 2x(b — a), where b is the radius of the outer conductor and 
a is the radius of the inner conductor. In the principal mode, the electric field lines are 
perpendicular to the inner and outer conductors, while the magnetic field lines are 
concentric circles about the axis of the inner conductor. When these conditions exist, 
the differential equations (1) apply. In the case of the tapered line, in order to satisfy 
rigorously the boundary conditions, other modes of propagation will have to exist. 
Fortunately, if the taper of the outer conductor is not too abrupt, it is possible to neglect 
the higher modes. Thus, this is another source of error between calculated and measured 
values. 

3. The inner conductor of the lossy tapered line is a glass tube coated with a metallic 
film. In our solution, we neglected to take into account the fact that the electric field 
penetrates through the film. As shown by Griemsmann,’’ for a uniform line, the effect 
of this field can be expressed as a correction factor in the characteristic impedance and 
the attenuating factor of the transmission line. He shows that the first order correction 
for the field inside the film leads to the following expressions for the characteristic 
impedance and the attenuation factor: 

Z. = R,{1 — iX — [(K? — 1)/4 log K — 1/2](8bX)*}'"”, 
(25) 
lr = i8{1 — iX — [(K? — 1)/4 log K — 1/2](gbX)*}'”, 


where 8 = 27/\, X = Rd/2xR, , K = a/b, and R = resistance per unit length, \ = 
wave length, Ry = characteristic impedance of the lossless uniform line, b = radius of 
the inner conductor, a = radius of the outer conductor. For the size of the transmission 
line and resistance used in this problem, the factor [(K° — 1)/4 log K — 1/2](@bX)’ 
can have a value of 0.05 at the lower frequencies but is negligible at the higher fre- 
quencies. The change in Z, and I could be interpreted physically as a change in the 
size of the outer conductor; therefore, the values in Table 1, indicating the linear and 
exponential line outer conductor radii, may be separated more than indicated. Since 
this phenomenon is more evident at lower frequencies, it may account for the further 
discrepancy between calculated and measured values at lower frequencies than at higher 
frequencies. 

In the design of a linear tapered line termination, a fairly flat and low VSWR can 
be obtained, if the following design constants are taken into account. 

1. k in the exponential of Eq. (4) has a value of 1.3/Z, where L is the length of the 

taper in meters. 
2. In Eq. (22) if 


L/xX = 0.62 and R,rA/2nZ.L = 0.2, 


*J. C. Slater, loc. cit., p. 162. 
J, W. K. Griemsmann, Design of a lossy co-arial transmission line, Doctoral Thesis, Brooklyn Poly. 


Inst., June, 1946. 
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where L = length of the taper in meters, \ = wave length in meters, R, = total resist- 
ance of the inner conductor = RL, Z, = characteristic impedance of the uniform line 
having an outer radius b and an inner radius a of the input of the tapered line, then 


T=>Z.. 


3. If 4r/kX = 6, we will find that the series for U and W can be approximated by 
exponentials. In addition, if R[exp (kL) — 1]/kZ, = 3.5, then 


0.96 <= (1 — U)/(1 — W) S 1.04. 


Since kL = 1.3, and exp (kL) = 3.667, we find that (1 U)/(1 — W) is nearly unity if 


Rises & 1.7. 


Thus it can be stated that to satisfy the lower frequency range with a low VSWR, 
the film length must be long (L = 0.6A), and the total resistance must be near the Z, 
for extremely low frequencies. In the upper range of frequencies, we should have a total 
resistance of 1.7 times the characteristic impedance Z, . Thus we must conclude for 
practical reasons, that it would seem more advisable to use the tapered line for fre- 
quencies above 1000 megacycles. Of course, if the restriction on the VSWR could be 
lessened to permit a VSWR of 1.2 or 1.3, then a larger range of frequencies could be 


covered by the linear tapered line. 


CAPACITY OF A PAIR OF INSULATED WIRES* 
By W. HOWARD WISE (Bell Telephone Laboratories) 


1. Introduction. This problem has been treated in elegant fashion by Craggs and 
Tranter.’’’’* Their first two papers employ a conformal transformation of the free space; 
the third paper works with charge distributions. They end up with an infinite deter- 
minant of value zero in which the unknown capacity appears in one element, and 
conclude that “satisfactory numerical approximations can be obtained by keeping only 
the first few rows and columns’’.* 

This note is written to remark that the end result of a straightforward attack with 
bi-polar coordinates is 


Y / f _s €0 a ~ .} 
C = & js. + ; log : di hy 3 (1) 


where C = capacity, ¢, = dielectric constant of air, or other material outside the 
jackets, e = dielectric constant of jacket material, s = interaxial separation, x = outer 





*Received Feb. 7, 1949. 

1J. W. Craggs and C. J. Tranter, The capacity of twin cable, Q. Appl. Math. 3, 268-272 (1945). 

2J. W. Craggs and C. J. Tranter, The capacity of twin cable—II, Q. Appl. Math. 3, 380-383 (1946.) 

3J. W. Craggs and C. J. Tranter, The capacity of two-dimensional systems of conductors and dielectrics 
with circular boundaries, Q. J. of Math. (Oxford) 17, 138-144 (1946). 
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radius of jackets, v = radius of conductors, and the k are to be obtained from the fol- 


lowing set of equations: 





. ee 
Fia. 1. 
mk, os 2k. ari 3k; i 4k, = 5k; — ORS ae L. 
—k, oe Mok,» a. 6k p_ 10k, was 15k; — 2a 1/2, 
—k, = Ak, 4 mks _ 20k, = 35k; — 2? a 1/3, 
(2) 
—k, = 5k — 15k; a mk, _— 70k; — 22? 1/4, 
—k, <2 6k. = 21k a 56k, + m;k; —— 2? oe 1/5, 
where 
— + € + (€ — €&)(v/x)”” (2)" _ (27 — 1)! 
" @ — & +(e + &)(v/x)” \x t\(r — 1)! 
= 1, — @r — 1I)l/rKsr — DI. 
The 7-th equation is 
— ~(r7+x«-—1)!, 
lke = = + Ye peer = Dir! - (3) 


If one keeps n rows and n columns in the above set of equations, the result is identi- 
cally the same as that obtained by keeping n + 1 rows and n + 1 columns in the in- 
finite determinant in Craggs and Tranter’s last paper. 

2. Discussion. The series in k converges slowly when the two insulating jackets are in 
contact. What is worse, if one stops with a small number of k, terms, say five, using (2) 
as written, the resulting >~°_, k, will be too small, not merely because there are only 
five terms in the series but also because each k is too small. The consequence of all this 
is that one should start out with the idea of using ten or twelve k, terms. Fortunately, 
this does not mean that we have to solve ten or twelve simultaneous equations by 
gradual elimination of the unknowns. Since the m are a rapidly increasing sequence of 
large numbers, the solution may be arrived at by successive approximations. 

If the wires were bare (x = v), the series in k would converge more rapidly; but in 


this case 
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(s/v)?? — (27 — 1)!/r'(r — 13, 


n 


>, k, = — log (1 — &,), 


n=1 


ky = {1 — [1 — (2v/s)"]"}/2, 


and (1) reduces to 
C = «/4 cosh™ (s/2v) 
(4) 


= € ‘4 log is Qu + [(s 2v)” _ 1)'7}. 
= 2, s/2v = x/v and 


When the two insulating jackets are in contact, s/x 


co kb (5) 


C = & 4 log 24° log x 
€E U n=1 


Ikpq 


0.6 0.7 0.8 0.9 t 





0.2 0.3 0.4 0.5 
u/Xx 


Fig. 2. 


The correction term >> k, in (5) is a function of ¢/¢) and v/x. Figure 2 is a plot of 7a, 
versus v/x for €/¢é. = 1, 2, 3, 4, 5 and . The curve for ¢/¢, = 1 was not obtained by 
= 6, then 


direct computation of the k but by noting that if « 
C = «,/4{logs/v — >> &,] 


€,/4 log {s/2uv + [(s/2v)” — 1]'”} 
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and so 


i k, = log 2 — log {1+ [1 — (v/x)*]"7}. 


The upper end of the curve for e/¢€, = 5 was fixed by noting that if v/x is close to unity 
then k, = ki/n, (m + Ik, = 1/0 — &), ky = {1 — [1 — 4/(m + 1)]'7}/2 and 
> k, = log 1/(1 — ki) = log 2/{1 + [1 — 4/(m, + 1)]'”}. 

With v/x = 0 it is easy to solve (3) by the method of successive substitutions, thus 
getting 


2n o 


k, = 9 e {1 + gl(po/p:)” — 1] + > g'[(p.-1/p.)” — (1/0, (6) 


‘= 


where 


g (e — €)/(e + €), 


Po = i. Dp, = te p; P2 = oe 2p’, tee P. = Pi-i — p D.-2 ; (7) 


rs) 


> k, = —g logp, — p g° log (p.p.-2/pi-1). 


s=2 


With p = 1/2 we have 





dk, = — > g** log™ (8) 


n=2 


This formula was used to check the more laborious computations based on (2) at v/x = 0. 
Since 1, = (s/v)”’ if € = € , (7) with g = 1 and p = v/s gives , k, for bare wires. 
Thus 


, ] — 
= €& {log + log p, + ,é log pp.-r/v'-s} 
s=2 


a 
| 


II 


€,/4 i{log + log lim p,/p,- ‘ 


io 


(4) 


€o {log a log 5 {1 oo [1 ae aot} 


= ¢,/4 log {s/2v + [(s/2v)” — 1]’”’}. 
It is easy enough to write out a solution of (3) by the method of successive substitu- 
tions and so find that 


rs) 


= — ] n+K—- 
do kn = Set Lp Le 


n=1 n=1 “n «=1 


top Soret por gf a 


a «= 


de + Perr ‘yo has + ees 


‘. h=1 dh w=l uh 
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but, aside from the simple case v/x = 0 discussed above, (9) seems to be of no help. 
It is equally easy to write out the Fredholm type of solution of (8), 


] ll 
Lk = +5 De Au, (10) 


k=l] 
where 


K., & | 


iK 


a a 


Ke KK | 


K. Kk. £= 


tm 


s.=K,-> ‘l+#eE5|K. & Kn} --, 


Knew Kus Kam 


and 


but no use has been found for it. 
3. Numerical example illustrating convergence. 


s = 2x, v/x = A, e/eé, = 4. 


If we assume that 1, 2, 3, 4 and 10 k& terms are enough, we compute in turn, k, = 
2097, >-2., k, = .2544, 503. k, = .2692, >>*_, k, = .2751 and > °°, k, = .2807. The 
sum of the first four of the k, in p ‘, k, is .2797. The late terms in the sum are more 
important indirectly in computing the earlier terms than in the sum itself. 


THE TEMPERATURE IN AN ACCRETING MEDIUM WITH 
HEAT GENERATION * 


By A. E. BENFIELD (Cruft Laboratory, Harvard University) 


The thermal problem of this note was solved with the hope of using it to try to test 
the theory that the earth was formed by accretion on the dust cloud hypothesis;’ but 
there are many uncertain and unknown physical factors involved and, on reflection, it 
seems that the contemplated thermal considerations are unable at present to help in 
drawing conclusions. However, it is hoped that the following mathematical solution 
may be of interest and aid to others having related problems involving less uncertain 
physical conditions. 

As the spherical case presents some difficulties, we shall merely consider here a 





*Received March 25, 1949. 
F. L. Whipple, Scientific American 178, 34 (1948). 


1See, for instance, 
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homogeneous semi-infinite medium moving with constant velocity v in the positive 
x-direction, as shown in Fig. 1. Material is supplied both at a constant rate and at a 











a Vv _ 
4 ff 
vr, 
tA 
y ie 
3s ¥Y 
Z fo 
‘Al 
Z , woo tet, 
2v¥ 
Z 
F Z a/V 
t=0 
of 
Z 
4, 
Z 
Vg ff I 1 j 
I 
o 2 3 4X 
Fic. 1. Curves of 7'(z, t) as given by (5) fort = 0,t,, and& = 4t;, when a, v, K and ¢, all equal unity. 


The plane x = 0, which is held at a fixed temperature 7 , represents the free surface of the accreting 
semi-infinite medium at t > 0. The dashed line of slope a/v represents the asymptotic value of the 
temperature for large values of t. 


steady temperature 7, to the surface x = 0 of the medium, in such a way that the 
medium always exists continuously in the region 0 < « <@, while for z < 0 there is a 
vacuum. In other words, x is measured to the right from the free surface of an accret- 
ing semi-infinite homogeneous medium. f 
For t < 0 and 0 < x <o the temperature is everywhere a constant, JT = JT, . At 
time t = 0, however, the generation of heat begins throughout the semi-infinite medium 
at the constant rate A cals/em*-sec and this heat generation continues indefinitely. 
For t > 0 the surface x = 0 continues to be held at the constant temperature T = T,, 
and the material accreting at the plane x = O contains the same uniform distribution 
of heat sources, A cals/em*-sec, as the rest of the material in the region 0 < x <@. 
The problem, then, is to find 7(2, t) for this situation when t > 0 and 0 < x <o&. 
The appropriate differential equation is 
fa KES - E+, (1) 
ot Ox Ox cp 


where c and p are, respectively, the uniform specific heat and the uniform density of 





+The accreting material arrives at the plane x = 0 with zero velocity, so that it carries no kinetic 


energy, and hence no heat is evolved on its arrival. 








438 NOTES [Vol. VII, No. 4 


the homogeneous medium, and K is its diffusivity (defined as the thermal conductivity 
divided by cp). Accordingly, K, v, A, ¢ and p are all constants in this problem. 7 
Using the method of the Laplace transform we now introduce the quantity 7’ = 
L{T} = fo eT (a, t) dt, where p is a constant. Defining a = A/cp we may follow the 
usual procedure of the Laplacian method and rewrite Eq. (1) as 
Sy 4 dT — : 
K v — pl +T, + a/p = 0, (2) 
dx” dx 
which is a second order differential equation with constant coefficients. The solution 


of Eq. (2) may be written as 


T = + 3+ ce" A, exp {a[(v?/4K’) + p ‘Ky'”"} 
(3) 
+ A, exp {—a2[(v’/4K’) + p/K]'”}], 
where A, and A, are functions of p, found by inserting the boundary conditions into 
Eq. (3). Proceeding now to do this, we may say that (i) when « = 0, T = JT, and T = 


T./p, from “pa : follows that a/p = —(A, + A,); and that (ii) ast @, T—> T+ 
at, so that T — »/p) + a/p’, from which we have A, = 0. 
Combining a ‘se values of A, and A, with Eq. (3) we find that 
— To a / 2) v2/2K { T(os2 72 {yz 1/2 
f= ; + D — (a/p')e exp {—a[(v'/4K°) + p/K]"*}. (4) 


The inverse transforms of the first two terms are 7, and at, respectively, and, as shown 
in a recent paper, 


exp {—\ (yp? /4K’) + p/K)'"} _ rf! E =/2K ' (z+ es) 
L\5, (x + vt) erfc Ky” 


~vs/QK;, ow ae, = oF 
a (a vt) erfe (2; =) |} 


ev 


erfe y = 1 — erf y = 1 — 2)” | exp (—u’) du. 


/0 


P 


where 


Our solution may now be written by adding the inverse transforms of Eq. (4); it is 


ie _ om _ @ | o2/K;_, eer zt) - a (Za%) | = 
T(x, t) = T, + at oy E (a + vt) erfe (22 We (a vt) erfe KD? (5) 


Curves of this expression for three values of the time are shown in Fig. 1. Equation (5) 
is seen to reduce properly for the special cases of x = 0, a = 0 and ¢ = O. Furthermore, 
it satisfies the differential equation (1), and it reduces as v — O to the appropriate ex- 
pression for a a sti utionary medium.° 





. Benfield, Q. Appl. Math. 6, 439 (1949), in which Eq. (11) is equivalent to the transform given 


here. 
3Cf. H. S. Carslaw and J. C. Jaeger, Conduction of heat in solids, Clarendon Press, Oxford, 1947, p. 60, 


where Eq. (2) is for a similar but not identical case with v = 0. 
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A point of some interest is the fact that, as time goes on, the slope dT (z, t)/dx of 
the curve of Eq. (5) does not increase without limit, but instead approaches an asymp- 
totic value a/v, shown in Fig. 1. This asymptotic behavior may be seen by differentiating 
Eq. (5) with respect to z, and then inserting the conditions ¢ > x/v and t > K/v’. 

It is a pleasure to thank Professor P. Le Corbeiller, Professor F. Birch and Professor 
I’. Whipple for helpful discussions. 


DEGENERATE TWO-DIMENSIONAL NON-STEADY IRROTATIONAL FLOWS 
OF A COMPRESSIBLE GAS* 


By N. COBURN (University of Michigan) 


1. Introduction. A class of non-steady, two-dimensional, irrotational, compressible 
flows which are very similar to steady, two-dimensional, irrotational, compressible 
flows will be studied. In order to do this, we introduce the well-known potential equation 
and Bernoulli relation for general non-steady flows. Our degenerate flows are defined 
by requiring that two families of cylindrical characteristic surfaces (with generators 
parallel to the time axis) exist in space-time. These flows have the following properties: 
(1) the wave fronts are stationary; (2) each of the velocity components and the speed 
of sound depends upon a single function of time multiplied by appropriate functions, 
which we shall call ‘‘reduced”’ velocities, of the space variables; (3) the single function 
of time is such that the motion decays as time increases. A canonical characteristic 
system, consisting of five equations with five dependent variables (the reduced velocities 
and the rectangular coordinates of the plane) and two independent variables, is ob- 
tained. It is shown that simple waves do not exist. Finally, it is shown that a degenerate 
non-steady flow, whose stream lines are logarithmic spirals, exists. 

2. The system of flow equations and the potential equation. Let x*(A = 1, 2)de- 
note a rectangular Euclidean coordinate system in the physical plane, and let ¢ denote 
the time variable. If v (A = 1, 2) denotes the components of the velocity vector in the 
x-coordinate system, and p and c denote the density and local speed of sound, re- 
spectively, then the equations of motion and the equation of continuity may be written 


as 





Ov, OV c” ap 

+ v* — ——=}0 (2.1) 
at Ox" P ax” F 
Op av” x Op 
g PE =@. 2.2 
nT Pa (2.2) 


In the above equations, the contravariant and covariant components of a vector are 
equal since the coordinate system is Euclidean orthogonal. However, we have introduced 
the notation of tensor analysis in order to use the summation convention. 





*Received Jan. 31, 1949. This work was done under a University of Michigan contract with the 


U.S. Army Air Forces. 
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For irrotational motion, a velocity potential ¢(¢, x”) exists such that 


Op ‘ 
n=, (2.3) 
Ox 


Further, we shall define the time component of the generalized velocity to be 


Ob 
—— 2. 
Vo at (2.4) 


The integrability conditions for (2.3) and (2.4) are the irrotationality conditions 


Ov, -_ Ov, Ov, = Oo (2.5) 
Ox" dx” dat da” 


By substituting (2.5) into (2.1), we obtain the Bernoulli equation 

2+ 9° + 2P = 0, (2.6) 
where g’ = vv», , and ¢ has been chosen so that the right-hand side of (2.6) vanishes. 
For adiabatic flows P is defined by 


p Cc Cc 

P= [d= —, y= 14. (2.7) 
J p q=% 

By differentiating (2.6) with respect to t, x‘ (where P is defined by 2.7) and substituting 


into (2.2), (2.1), we obtain the desired potential equation 


Ov, , Ov Ou ; 
a a? +a" — oO, (2.8) 

Ot Ot dx" 

Nu Nu 2 hu Awe P 

where a” = v'v" — cg” and g”™ is the metric tensor. 

3. Degenerate flows. Let w(t, x‘) = constant denote the equation of a family of 
o' characteristic surfaces. These loci satisfy the first order equation’ 
Ow \" , Ow Ow Ow Ow , 
(2 ) + 2° —— + a” s—— = @, (3.1) 
dl dx” Ot Ox” Ax" 


If a flow possesses two families of characteristic surfaces which are cylinders with gene- 
rators parallel to the time axis (that 7s, w is not a function of t), then the flow will be called 
degenerate. 

For degenerate flows, (3.1) reduces to 

1. OW Ow 
mu Ow Ow _ 29 
a ‘= 0. (3.2) 


OX 


a 


Ox 
From this fact, it is easily verified that these flows have the following properties: 
(1) they are always supersonic; 
(2) the normal vectors to the two families of characteristic surfaces and the corre- 
sponding bicharacteristic vectors are independent of time; 
the direction of the velocity vector and the Mach number are independent of 


(3) 


time; and 


IR. Courant and D. Hilbert Methoden der mathematischen Physik, vol. I1, Julius Springer, Berlin, 


1937, p. 375. 
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(4) if »*, ’n* denote the unit normal vectors to the two families of characteristic 
. x ° . ° ‘ot 
surfaces and ?*, ’t’ denote the unit vectors of the corresponding bicharacteristic 


curves, then 
y -_ cn* 4. (q° —_ fy". (3.3) 
v=c'n — (¢ — cy)". (3.4) 


From property (2), we see that the wave fronts are stationary. In view of property 


a3) 


(3), we may write 
y = f(t’, c = f(ie, (3.5) 


where the barred quantities are the “reduced”’ velocities, ete., and depend only on the 


variables z*. From (2.6), (2.7), it follows that 
vo = [f(D] % . (3.6) 


Substituting (3.5), (3.6) into (2.5), (2.8), we find that f(é) must be 
: I - 
JO = -y TF (3.7) 


where k, b are constants. Further (2.5), (2.8), (2.6), reduce to 


Ov, Ov OVo - ' 
potest — —_ = kp (3.8 
Ox" da" da i oe 
Ov) gs _9 
a“ * + Qk, + q) = 0, (3.9) 
Oa” 
uy ‘ 2¢” 
BW +g t+ ere | = 0. (3.10) 


In our future work, we shall drop the bars over the ‘‘reduced velocities” etc, but shall 
limit our considerations to the system (3.8) through (3.10). 

4, The characteristic system for degenerate non-steady flows. In view of (3.8), 
Eq. (3.9) is a quasi-linear partial differential equation of the second order. The char- 
acteristic system for such an equation can be obtained immediately from H. Lewy’s 
work.’ Let a = constant, 8 = constant denote the traces of the two families of cylindrical 
characteristic surfaces on the x’-plane. In addition, in order to compare our results 
with those of the steady case, we introduce the notation 


2 1 2 
z= 2 y=2 “u=y=0, v=%,=0. (4.1) 


An application of Lewy’s method (and a little algebra) furnishes the characteristic 


system 


2R. Courant and D. Hilbert, loc. cit., pp. 326-332. 
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: Sade a) > are Ox : 
(u(q’? — c’)'”? + ov) > — (v(¢’ — c’)'” — cu) a = 0, (4.2) 
(u(q — c)'”? — cw) : © ew: (gq — c’)'? + cu) : Ems 0, (4.3) 

Oa 0a 


0 2k Fi 7) F w 
du 4 kg + q) 92 = 0, (4.4) 


ai 1/2 Ov 2 2\1/2 
‘ — + 7) a —_— — (PY ‘ a 5 
(v(q c’) + cu 98 + (u(q c) cv) ap ' ug —c)* +c oB 


as 1/2 Ov 2 2\1/2 Ou 2kq'(v + q) ox 
(a — cy”? — cu) — ila =e) ”) = 5 Ft = 0, (4.5 
(9 . mm 0a + (ug 7 ¥ @) 0a T ug — sos — cv 0a ») 
0 ’o 03 g ) 
> — ku : Sas kv oy = (). (4.6) 
0a Oa 0a 


The algebraic relation (3.10) must also be added to the system (4.2) through (4.6). 

Somewhat different forms of (4.4), (4.5) can be obtained by solving (4.3), (4.6) for 

dx/da, and (4.2) and the equation in 8 corresponding to (4.6) for 07/08. We find that 
dx _ (u(g’ — c’)'” — cv) A dx (u(g’ — c’)'” + cv) du 


Oa k¢(q —c)'" da 0p kq(q —c)'” op a 


When (4.7) are substituted into (4.4), (4.5), we obtain equations involving derivatives 
of u, v, v . Further, if we introduce 6, the slope of the velocity vector 
u = q cos 8, v = qsin 8, (4.8) 


then (4.4), (4.5) become 


i n1/2 O@ 2 00 [Avo + gq’) Io 
(a —c - q = a = 0 (4.9 
qd ee 0B veg 0p = (¢ —c)'” 0B : } 
2 2\1/2 0q 2 06 (vo + q’) Oo 
(a —c —+~ — cq - ; ss. (= |G. t.10) 
qq c ) Bax CC Ja + ie = c’)'? rw ( 


Simple waves, (that is, a family of straight line bicharacteristics, say 8 = constant, 
along which u, v are constant) do not exist for degenerate non-steady flows. For if we 
require u, v to be constant along a curve 8 = constant and the slope of (4.3) to be in- 
dependent of a, then a differentiation shows that c must be constant along 8 = constant. 
From (3.10), it follows that vp is also constant along 8 = constant. Use of (4.5), shows 
that ¢g = 0 or »% +g = 0 along curves of this family. Finally, equating the slopes ob- 
tained from (4.6) with v) constant and (4.3), we see that g = 0 or g = c. Hence g = 0 
is the only common solution of our equations. 


The case k = 0, vo = constant, furnishes steady two-dimensional flow. Another case 
of interest is when vo + q° is zero throughout the flow. From (3.10), it follows that such 
a flow can exist only at the Mach number, M = 5'”*. The Eqs. (4.9), (4.10) may be 


integrated immediately. We obtain 


@ tan 


q = a "a along a = constant, (4.11) 


6 an w@ ° ‘ 
_ along 8 = constant, (4.12) 


q = be 

















1950] C. C. LIN 443 


where w is the Mach angle for M = 5'” and a, b are arbitrary functions of a, 8, re- 
spectively. The curves (4.11), (4.12) are logarithmic spirals in the hodograph plane. 
For the case v% + g = 0, the (z,y)-plane map of the bicharacteristics a = constant 


(or 8 = constant), is orthogonal to the (u,v)-plane map of 8 = constant (or a = con- 
stant). Hence, the (x,y)-plane map of (4.11), (4.12) consists of logarithmic spirals. Specific 
formulas can be obtained easily. Since the bicharacteristics are inclined to the stream 
lines at the constant angle w, the stream lines are also logarithmic spirals. Evidently, a 
similar method can be used to investigate those degenerate non-steady flows for which 
Vo f q),S(7) < @. 


NOTE ON THE CHARACTERISTICS IN UNSTEADY 
ONE-DIMENSIONAL FLOWS WITH HEAT ADDITION* 


By C. C. LIN (Massachusetts Institute of Technology) 


1. In a very interesting paper, Kahane and Lees’ studied the problem of one- 
dimensional wave propagation in a gas when heat is being added. They used the method 
of numerical integration, by taking finite differences along the characteristics. However, 
since the form of the characteristic equations contained more than two dependent 
variables, they were led to use some rather artificial approximations besides those in- 
volved in taking finite differences. As the type of work done by Kahane and Lees will 
probably be continued by people interested in jet propulsion, it seems desirable that 
the most convenient form of the characteristic equations be derived and a simple pro- 
cedure of numerical integration be developed. It is the purpose of this note to show 
that with a proper choice of the dependent variables, the characteristic equations are 
much simpler, and the numerical integration can be carried out in a straightforward 
manner. 

2. The fundamental equations for one-dimensional unsteady flow are 


“ue + wu, + Se = 0, (1) 
p 
p: + up. + pu, + puA ‘A, = 0, (2) 
s, + us, = q*, (3) 
p = constant X pe’, (4) 


where u is the velocity in the direction of the z-axis, A is the cross-sectional area of the 
tube, ¢ is the time, and p, p, s, y are the familiar symbols for the pressure, the density, 
the entropy (divided by the specific heat at constant volume c,) and the ratio of specific 

*Received March 21, 1949. 

14. Kahane and L. Lees, J. Aero. Sci. 15, 665-670 (1948). 

?For further discussions of this problem, see a note by William Swartz, which is to appear in J. Aero. 
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heats, g* is the rate of production of entropy, i.e., it is the heat generated per unit mass 
Y J? b 
divided by the product of absolute temperature 7 and the specific heat at constant 


volume ¢, . 

The Eqs. (1)-(4) form a system of three differential equations in three variables, 
provided one of the four variables p, p, s, u is eliminated by the use of (4). One can, 
e.g., eliminate the derivatives of p from (2) to get 


Di + up. + ypu, = p(q* — yuA~"A,). (2A) 


The advantage in this choice lies in the fact that (1) and (2A) contain only the de- 
rivatives of the two variables p and u. 
3. The standard method of finding characteristics can now be applied. If we write 


p, dt + p, dx = dp, 
u, dt + u,dx = du, 
s, dt + s,dx = ds, 
and combine these equations with (1), (2A), and (3), we have six equations for the six 


partial derivatives. It can easily be seen that the coefficients and the right-hand sides 
of these equations can be arranged into the following matrix: 


(1 uw 0 p? 0 0 O 
10 yp 1 u 0 0 p(q* — uyA"*A,) 
109 0 O O Ll ui qf 


idt dx 0 0 0 0 du 


10 O dt dx 0 O dp 





0 0 0 0 dt dx ds ) 


The coefficient determinant breaks down into two distinct factors 


dx — udt = 0, (dx — udt)’? — a’ dt’ = 0. 


As was expected, these are the familiar characteristic base curves. The corresponding 


compatibility equations can easily be obtained: 


ds 


lx t- 
7 q* along — = u, (9) 


dt 


dp du q* : dx s 
ot +—= (< —uA A.) dt along =- =~uta. (6) 
yp a Y dt 
One may note in passing that these specialize immediately into those in the case of 
uniform entropy, provided the variables p, u, s are also used in the latter case. In that 
case, with A, = 0, it is possible to integrate the compatibility equations (6). 
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4. In the process of numerical integration, to reach a point P; from a knowledge 
of the solution along the curve P,P, (Fig. 1) one may use the set of equations (6) to 
Uz, and ps , since g* and A are supposed to be known functions. The 


Ps 


calculate z; , ts , 


Fia, 1. 


path line P,P; can then be drawn in accordance with the second equation of (5); the 
first equation of (5) will give the value of s; . One can then calculate p, from (4) and 
obtain all of the dynamic and thermodynamic variables at P; . Iteration processes can 


be carried out in the usual manner. 


A NEW SUPERPOSITION PRINCIPLE FOR 
STEADY GAS FLOWS* 


By R. C. PRIM** (Naval Ordnance Laboratory) 


This paper is concerned with steady flows in the absence of extraneous fields of 
force of a frictionless, thermally-nonconducting gas having a product equation of state, 
i.e., an equation of state connecting density, pressure, and specific entropy in the form 


p= P(p)S s). 
H. Poritsky [1]? has discussed the construction of steady, spatial gas flow solutions 


from steady plane flow solutions by the superposition of a uniform velocity field normal 
to the given plane flow field. In particular, he pointed out that if 

V, = iu(a, y) + jo(a, y) (1) 
is a plane velocity field (referred to ordinary rectangular Cartesian coordinates x, y, z 
with unit vectors i, j, k) satisfying the equations of steady-state gas dynamics, then 
the spatial velocity field 


V=V,+Kk),, (2) 


where V,, is a constant, also satisfies those equations. 

*Received March 21, 1949. 

**Now with Bell Telephone Laboratories. 

tNumbers in square brackets refer to the bibliography at the end of the paper. 
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The validity of this superposition principle follows at once from (Newtonian) rela- 
tivity considerations. In reference to an observer having a uniform velocity —kV, 
with respect to the flow field (1), this field has the velocity vector (2). It should be 
noted that it is essential that the given field (1) not depend on z, otherwise the flow 
field with respect to the moving observer will be unsteady. 

The basic superposition method of deriving spatial flow fields from a given plane 
flow solution will now be given a new order of power by making use of the Substitution 
Principle for steady gas flow. 

This Substitution Principle was first established by Munk and Prim [2] for the 
vase of gases having constant specific heats, and was extended by the present author [3] 
to the broader class of gases having a product equation of state, p = P(p)S(s). It can 
be stated briefly as follows. If V, p and p are, respectively, the velocity, density, and pressure 
of a possible flow of a given gas, then mV, p/m and p are the corresponding quantities of 
another possible flow of the same gas, provided only that V-grad m = 0, that is, that m is 
constant along each individual streamline. 

The members of a family of flows related by this Substitution Principle clearly share 
the same streamline pattern and pressure field and, as may easily be verified, they have 
also a common reduced velocity field W. (W is defined as V/a, where a is the ultimate 
velocity magnitude attainable on a given streamline by expansion to zero pressure, 
a’ = V* + 2h, with h the specific enthalpy.) In terms of a given reduced velocity field, 
the Substitution Principle manifests itself in the arbitrary assignability of the ultimate 
velocity a for each streamline. Hence, it is the reduced velocity field W rather than the 
actual velocity field V which plays the basic role in the theory of rotational gas flows. 
We, therefore, focus our attention on the problem of obtaining spatial reduced velocity 
fields satisfying the equations of gas dynamics (see [4]) from plane fields that do so. 

Specifically, we suppose given a plane flow field with reduced velocity W, and inquire 
what spatial reduced velocity fields W can be obtained by a combined application of 
the Substitution Principle and the Newtonian relativity considerations discussed above. 
We write W in the form 


W = aW, + kW. (3) 


Denoting the ultimate velocity functions for the plane flow field and the spatial flow 
field by a, and a, respectively, we have: 


V, = a,W,, (4) 
V = aW, (5) 

and 
\. = ay... (6) 


The ultimate velocity functions a, and a are defined by the following relations: 
a, = V,+ 2h, (7) 
a’ = V* + Qh. (8) 


Now the specific enthalpy of the plane flow is not affected by the superposition of a 
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constant V, . Hence h, = h. Furthermore, from (2), V? = Vi + V: . Equations (7) 
and (8), therefore, yield the following relation between a, , a, and V, : 


a, =a’ — Vi, 
whence, making use of (6), 
a, _ a 72\1/2 
7 = (1 — W2) (9) 
and 
Wi a ees (10) 
." W+aQ™ 
From (2) and (5) we obtain 
. 
W= = +k Vn 
a a 
(11) 
7. 
= Vee yy Ve 
a, a a 
Substitution of (4), (6), and (9) into (11) yields 
W = (1 — W,)'?W, + kW, . (12) 


Equation (12) is a formula for constructing spatial reduced velocity fields from a 
given plane reduced velocity field and a function W, . The nature of this function W, 
is clarified by (10) and the Substitution Principle. V, can be any (positive or negative) 
constant, while, by the substitution principle, the reduced velocity function a, can be 
assigned any (positive) values constant along each individual streamline without 
affecting W, . Therefore, letting y,(zx, y) be a streamfunction of the given plane field, 
it is only required that in the spatial field the surfaces of constant ¥, coincide with the 
surfaces of constant W,, and that W,, be of one sign throughout. That is, 


W, = Walvo(z, y)I; (13) 


where W, is either a non-negative or a non-positive function everywhere. In vector 
form, the restriction (13) can be written 


W, - grad W, = k - grad W, = W - grad W, = 0, (14) 
or 


OWn _ 


ap = 0 (15) 


W, - grad W, = 
Of course, physical meaningfulness requires also that Wi < 1. 

In the formulas (12) and (13) or (12) and (15) we have a means for generating a 
vast variety of spatial reduced velocity fields of possible gas flows from the reduced 
velocity field of any given plane flow. The power of this generating method stems, of 
course, from the essentially arbitrary nature of the function W,|y,]. 

One question that naturally presents itself is whether the generating formula (12) 
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is rotation-preserving; that is, does curl W, ~ 0 imply curl W + 0? To answer this 
question we compute from (12): 


W 


curl W = (1 — W;)'” curl W, — ies Tay sgrad W, X W, + grad W, X k. = (16) 


From (16) we find that 
k X curl W = grad W,, 


from which it follows that curl W = 0 implies grad W,, = 0. But, together, curl W = 0 
and grad W, = 0 imply (from (16)) that cuyl W, = 0 (for W; < 1). Therefore, curl 
W is zero if, and only if, both curl W, and grad JW, are zero. Hence rotational W fields 
san be derived from irrotational W, fields, but the converse is not possible. 

Another class of spatial flow fields of particular interest are the “generalized Belt- 
rami flows” investigated by Nemenyi and Prim [5, 6]. These are the flow fields for 


which 

W X curl W = 0 (17) 
throughout. For plane and axially-symmetric flow fields this condition is evidently 
equivalent to the condition, curl W = 0; however, for more general spatial flows, the 
class of flows satisfying (17) is much larger than that satisfying curl W = 0. Aside 


from their interesting geometric and kinematic properties, these generalized Beltrami 
flows have a special physical significance, shown by the relation 


grad log H(p ‘= 2 W seu : 


(18) 
where H(p) = f¢ [1/P(p)] dp and p, denotes the stagnation pressure (pressure on a 
given streamline for zero velocity). From (18) it is seen that (for W* < 1) the class of 
generalized Beltrami flows is identical with the class of flows having a uniform stagnation 
pressure. The given general form of (18) is due to the present writer. The particular 
form valid for gases having constant specific heats (and hence a constant adiabatic 


exponent y) was published earlier by B. Hicks and his colleagues [6]: 


y= 4 W X curl W 
grad log po = . 


= aA (19) 
“y 


(For the case of gases having constant specific heats, the function P(p) can be taken 
as p’’” whence H(p) = {3 p''’” dp = y/(vy — Ip "””.) 
We shall now make use of this property of generalized Beltrami flows to establish 
a method of generating them through the Superposition Principle of the present paper. 
In addition, we shall need the following two relations: 
H(p) = H(p.)(1 — W), (20) 
and 
(1 — W’) = (1 — W,)(1 — W). (21) 
The latter follows at once from (12). The former is the general form (ef. [3]) for gases 


with a product equation of state of the familiar formula 


p= p(l — W’)""" (22) 
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assuming constant specific heats. Now the application of either the Substitution Principle 
or the Newtonian Relativity Principle does not affect the pressure field of the given 
plane flow field. Therefore, denoting the stagnation pressure of the given plane field 
by Pop , We have from (20) that 


H(p) = H(po)(1 — W*) = H(por)(1 — W2); 
making use of (21), we obtain 
H(p.)(1 — Wi) = H(p,). (23) 


(It should be remarked here that po, , like W,, , is a function of ¥, .) Now H(po,) is a 
known function of x and y once the plane field is given. It can be computed from a given 
W,, by integration from (18). Therefore, in order to obtain from (12) a spatial field 
which is a generalized Beltrami field, it is only necessary to choose the function W,, so 
that in (23) H(p») is a constant. Since 0 < W < 1, the constant must be chosen so 
that throughout the region of the given plane field considered 


H(po,) < H(po). (24) 


This implies that any curves in the plane flow field on which H(po,) = © must be ex- 
cluded from the region considered. 

We thus have a simple method of constructing a one-parameter family of generalized 
seltrami flow fields from any given plane flow field: H(p.) is simply assigned any con- 
stant value satisfying (24), and the function W,, for use in (12) is then computed from 


a Pe Hon) |" , 
y= E H(pa) . (25) 


It should be noted that if curl W, = 0, (18) implies that H(po,) is a constant, so that 
(by (25)) setting H(po) equal to a constant forces W, to be constant also; hence by (16) 
curl W = 0. Therefore, only rotational plane fields generate spatial fields for which 
W X curl W = 0 while curl W = 0. 

This powerful method of generating generalized Beltrami flows produces, as quite 
special cases, all previously known examples of such fields. 

As an example of an application of this method producing new Beltrami fields we 
consider the following one-parameter family of plane rotational flows, discovered by 
the present writer [7], for gases having constant specific heats: 


. = 6 
W, = 1, sin — + 6, cos — + Z,-0 (26) 

NUo Ao 
where r, , 6, , Z, are unit vectors in an ordinary cylindrical coordinate system (r, 0, Z), 
h = (y + 1)/(y — 1) and y is an arbitrary constant such that 0 < v5 < 1. For this 

set of fields we compute (from (19)) 

6 (A+1) (049—-1/A)/ (1-070) is 
Pop = (; cos" o) (27) 


and obtain the two-parameter (v) , A) family of generalized Beltrami flows: 
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‘ 6 gyre 
W r,A sin (, cos. ) 


NUo \ NVo 


150 


4 9 \(ote- I= 00) 
+ 6,Av, cos (; cos’ ) (28) 
Avo 


+ z,| 1 _ A'(r cos’ 2) ; | a 


where A is a constant such that 


6 (v29—1/A)/(1—0%0) 
G<A2:< (; cos ) (29) 
No 


throughout the flow region considered. The lines in the 7, @ plane on which r cos* (@/Av) 
is equal to a constant (for fixed vy) are the streamlines of the field (26). Therefore, for 


v5 1/X the innermost streamline in the region considered fixes the upper bound on 
A, and it is clear that the line 6 Avo /2 must be excluded. For vj > 1/d the point 
r =o js the critical one in bounding A; only a finite part of the plane can be taken as 


the region considered. 
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ON HEAT TRANSFER PROBLEMS IN VISCOUS FLOW* 
By G. F. CARRIER anp J. A. LEWIS (Brown University) 


1. Summary. Many problems of physical interest which are associated with the 
flow of a viscous fluid through a narrow channel require the determination of the tem- 
perature distribution throughout the field of flow. In general, such problems may be 
separated into one of three classifications. The first of these is characterized by the 
existence of a thermal boundary layer, the second by a temperature distribution in- 
dependent of the coordinate across the channel, and the third by an intermediate type 


*Received March 7, 1949. 
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of temperature distribution. In the present paper a simple method for the classification 
of such problems is given, and techniques for solving some of the problems which arise 
in practice are discussed. 

2. The general equations. We shall consider the two-dimensional flow of an incom- 
pressible fluid with constant viscosity and heat conductivity. The equations governing 


this flow may be written 


p(uu, + vu,) = —p, + pwAu, (1) 
p(uv, + w,) = —p, + wAr, (2) 
uz +v, = 0, (3) 


where u, v are the velocity components in the x, y directions, respectively, p is the pres- 
sure, p the density, » the viscosity, and A the Laplacian operator. 
To determine the temperature distribution, we also introduce the energy equation 


pc(uT, + vT,) = kAT + pl2uz + 2; + (u, + v,)7], - (4) 


where k is the heat conductivity coefficient, c the heat capacity, and T the absolute 
temperature. 

3. Channel flow. We shall consider the flow of a viscous, heat-conducting fluid 
through a channel bounded by smooth walls of small curvature, which may be moving 
or stationary. Then Eqs. (1), (2), (3), and (4) may be replaced by a simplified set of 
equations, which are obtained here by making certain physically plausible assumptions 
concerning the nature of the flow. We assume that u > v, p, > p, , and that the Reynolds 
number of the flow is such that the inertia terms in the equations of motion are very 
small. Under these assumptions the general equations can be replaced by the approximate 


° 1 
relations 


Uy = Pz = p’(2), (5) 
| u(x, y) dy = const., (6) 
kT,, — peuT, + wu, = 0, (7) 


where —?,(x) and ¢,(x) are the coordinates of the upper and lower boundaries. Typical 
boundary conditions’ are 


1A crude analysis will show that the term pcv7’, is of the same order as pcuT’, . However, since the 
channel is assumed to be bounded by walls of small curvature, the streamlines of the flow will practically 
coincide with the lines y = const. Therefore v will be taken to be zero to the order of approximation 
considered here. An equation having the general form of Eq. (7) may be obtained also by a suitable 
change of variable. For example, in the doctoral thesis of M. Finston (Brown Univ.) the details of this 
variable change for a specific problem are given and it is shown that the approximation used here is 
valid. 

2In general 7'(0, y) is given only roughly for the following reason: the wall effects near the entrance 
imply a modification of the upstream temperature distribution via conduction and thus, at any given 


section, 7’ cannot be specified a priori. 
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Pintee = Const., 

T (0, y) const. = T,, 

T(x, —t,) T(t, te) const. = T',, or 
‘= (T,)y-s, = const. 


It is evident from the form of Eq. (5) that the velocity distribution at any section 
(i.e., at fixed x) is parabolic. That this is the velocity profile one obtains in such problems 
is well-known. Thus we need consider only Eq. (7) with its associated boundary condi- 
tions and regard the velocity uw as known. 

It is convenient to make the definitions’ 


é = z/L, 

n = [2y + t(x) — 4(a)]/[t.(7) + t(a)], 
e = kL/pcU 6 , 

6 = (T — T5)pcts/pUoL, 


where L is the length of the passage, t) the minimum half-width of the channel, and 
U, the maximum velocity at the section of minimum area. Equation (7) then becomes 
(to our order of approximation) 
€6,, — [1 — )f(é) + UCL + n)/2 + UCL — n)/2]6, 
(8) 
+ (U, -—U,- inf)* i = Q, 


where f(£) is determined (in terms of p’(x), t.(x), and ¢,(2)) from Eqs. (5) and (6). 

4, The characteristic parameter «. Let us now examine Eq. (8) with regard to the 
quantity ¢, which depends upon the material properties of the fluid, the velocity of 
flow, and the geometry of the channel. 

If « > 1, Eq. (8) may be replaced by the approximate relation 


64, = 0. 


This implies that the temperature distribution across the channel, for «€ > 1, is inde- 
pendent of 7. There will, of course, be a thin transition region near the inlet where 
this is not true, but to the order of our approximation we can consider the temperature 
distribution independent of the coordinate across the channel in this case. 

When « is of order one, it is more difficult to compute the temperature distribution, 
since now all the terms in the differential equation are of the same order. However, it 
is still required that the temperature distribution should change smoothly as the fluid 
progresses down the channel. Thus a numerical procedure can be employed with reason- 
able efficiency to obtain the temperature distribution for any particular problem. An 
example of this type of problem is the flow of lubricating oil through a bearing.* The 

‘The definition of 7 is chosen to give 7 = +1 at the passage walls. 

‘For a particular lubrication problem, it was found that e 4. 








1950 G. F. CARRIER AND J. A. LEWIS 453 


usual assumption of constant temperature across the oil film is not exactly satisfied; 
however, probably little error will result in any practical case where this assumption 
is made 

When the heat conductivity of the fluid is small, the velocity of flow large, and the 
channel not too narrow, we will have e < 1. In this case we might expect that the first 
term of Eq. (8) could be neglected. However, if this is done, the solution found will 
never obey all the boundary conditions and in many cases will be singular at the walls. 
Hence, we must expect that the conductivity will play an important role near the walls. 
We shall anticipate that there will be a thermal boundary layer whose thickness is of 
order t,e” in which the temperature will vary rapidly from its value at the wall to its 
value outside this layer. In the latter region the temperature is almost entirely governed 
by the convective effects. From a crude dimensional analysis we might conclude that 
m ~ 1/2. We shall see later, however, that the precise value of m depends upon the 
boundary conditions imposed upon the velocity. 

The solution of Eq. (8), with suitable boundary conditions, will be found as the 
sum of two solutions—a solution of the homogeneous equation which takes special 
account of the rapid variation of temperature near the walls and hence is a boundary 
layer solution, and a particular solution of the original inhomogeneous equation. Let 
us denote the boundary layer solution by 6, and the particular solution by 6, . In general, 
the inhomogeneous solution will describe the convective contribution to the temperature 
field, whereas the homogeneous solution will describe the boundary layer effect. How- 
ever, when one or both of the passage walls is stationary, even the inhomogeneous 
solution must be subdivided into convective and conductive parts, or else the solution 
will be found to be composed of such complicated functions that computation becomes 
impractical. This will be seen in the example at the end of the paper. 

5. The boundary layer solution. To obtain the boundary layer solution, we consider 


the homogeneous equation 
€(8;)an — F(é, 0)(0:)¢ = 0, (9) 
where F'(£, ») is defined in Eq. (8). The flow near each wall is treated by itself, and the 
two boundary layer solutions thus found are combined to give a smooth fit at the middle 
of the channel. For example, suppose we consider the thermal boundary layer at the 
lower wall (i.e., at 7 = —1). We set 
6 = (n+ 1) a 


and choose m so that the convective and conductive terms arising from Eq. (9) under 
this substitution are of the same order in e. If only the terms of the lowest order are 
retained, / may be written in the simplified form 
F = U, + e§[2f(€) + (U2 — U:)/2). 

If the wall is stationary, U, = 0. If it is moving, U, will be the dominant term near 
the wall and the term in e” can be neglected. In the first case we choose m = 1/3 and 
in the second m = 1/2. Then Eq. (9) has the form (with p and g known functions of &) 

(O)ss — Sp(Ade = 9, (10a) 
if the wall is stationary, and 

(O:)r¢ — G(Ai)_e = 0, (10b) 
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if the wall is moving. These equations may be reduced to a simpler form by the trans- 


formations 
dé 
a - 
J plé) 


to the forms 
(Orr — C(O). = 0, (1la) 


(0;:)er — (0,). = O. (11b 
$5 


The function @, must now satisfy boundary conditions of the following type: 
é, = 0, fora 0 (at the inlet); 
6, — 0, asf — © (at the middle of the channel); 
6, = G(a), for ¢ = 0, 


where G is chosen so that 0) + @, satisfies the original boundary condition at the wall. 
The problem given at the end of the paper exemplifies the simple method of obtaining 
the boundary layer solutions by the application of the Laplace transform. 

The determination of a particular solution depends so much upon the form of Eq. 
(8) in a particular problem that few statements can be made, in general, as to its form. 
In many cases it may be found by reducing the partial differential equation thus found 
to an ordinary differential equation by an appropriate change of variable.’ 

6. The flow between two parallel walls. To illustrate the method of solution when 
e < 1, we shall compute the temperature distribution in the flew of a viscous fluid 
between two parallel, stationary walls. The temperature of the fluid at the inlet will 
be taken as 7’ 


(although, as previously noted, this boundary condition will not be 
satisfied exactly), and the temperature of the wall will be taken as fixed at a value 
T,, . In this case, the function F¥ has the form 


F(é, n) =l1- "; 
and Eq. (8) becomes” (rigorously in this example) 
€0,, — (1 — °)0 + 4° = 0. (12) 


By inspection, a particular solution of Eq. (12) is seen to be —n‘/3«. This solution, 
however, does not yield the proper behavior in the interior of the channel, where con- 
vection and dissipation are the dominant effects. The above solution takes into account 
conduction and dissipation, the effects of which are dominant only near the walls. We 
seek a solution which is of the order «€ '’* at the walls and which tends to zero at the 

See, for example, M. Finston, Doctoral thesis, Brown Univ. (forthcoming). 

‘We note here that the exact solution of this problem could be expressed as an infinite series of 
certain integrals of Whittaker functions, but this fact is obviously of academic interest only. 
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middle of the channel. To find such a solution, we split the particular solution @, into 


two parts, 6* and 6**, satisfying the equations 
e0*, — (1 — 9) = —4n7(1 — 7°”), (13a) 
e9** — (1 — 7 )0F* = —4n’”’, (13b) 


where v is some large number. If now Eq. (13a) is replaced by the approximate relation 
(1 — 9°)0f = 4n°(1 — 9 
6* will be given by the expression 


6* = 4tn —  o 


which vanishes for 7 = 0 and approaches the value 4 v é for large v as 7 approaches +1. 
If we had merely neglected the first term in Eq. (12) and then computed the corre- 
sponding solution, we would have obtained a solution which becomes singular at 7 = +1. 
This seemingly artificial device used above avoids this difficulty. To obtain a solution 
ol iq. (13b) we set 

(eo'*1+7), for-1S <0, 

le "(1 — n), forO0 Sn Sl. 
If now we retain only the terms of lowest order in e, Eq. (13b) becomes 

o¢* — 2rok* = —4€ °° (1 — € 5)”. (14) 


1 


* . . . a" -—1/3 . a : ’ 
Since e < 1, a convenient choice for v is 2v = e€ ’°. It is known that 


s\ 2” 
. é _ at 
lim (1 - ra =e, 


and thus, since we have assumed e very small (i.e., 2v very large), we may approximate 
Kq. (14) by 
** 9+-o** — —a7e_~f - 
An obvious solution to this equation is 
ar* _ aig ae * 
This solution tends to zero for large ¢(i.e., at the middle of the channel) and assumes 
the value —4e '’* at the wall. In terms of 7 then the complete particular solution 4 


has the form 


2 (1 — ” *) —2y + —2v(l1—7n 
6, = 4&n° —s x — 8r{e" ee ee 3 
a= 
where 2v = e '’*. This function takes on the approximate value 4r( — 2) at the walls 
and exhibits the proper behavior at 7 = 0. Furthermore, it practically vanishes all 


along the inlet section except near the walls. 
To complete the determination of the temperature distribution, we must now obtain 
a boundary layer solution, 6, , which satisfies 


(Oi)s¢ — 25(01:), = 0, (16) 
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and the boundary conditions 
0,(é, 0) = 4r(2 — 4) + 8, 
6,(&, o) = 0, 
6,(0, 5) ~ 0, 
where 8@,, is the dimensionless wall temperature. 
To obtain this solution, we set 
6(é, 9) = [ e *'g(t, ¢) dt, 
i.e., we assume that 6, is the Laplace transform of some function g. Then g must satisfy 
the differential equation 
9sr(t, $) + 2gtgt, 4) = 0, 
which has the solution 
g = A(t)E(9t/2¢*)'°1,3(2(2t°/9)'”) + Bt) (9t/2¢°) 7° J _13(2(2t¢°/9)'””). 
If A and B are taken as constants, the function 6, has the form,’ 
6,(é, ¢) = Age *” exp (—2¢°/9) + BE?” exp (—2¢°/9E). 


The first term may be written as, 


ae~4/3 f . < i a 3/ 
Até exp (—2¢°/9£) = 3.4 9 | exp (—2a°/9) da. 
OF Jeyerse 


Equation (16) implies that if 6, is a solution then (@,); is also a solution and vice versa. 


Hence P(é, £) is a solution, where 


<n 
(Ee F&F) = - ~ 4 — 9,,3 /Q) 
PE, §) T(2/3) | = exp (—2a°/9) da. 
The quantity (2/3) is the Gamma function with argument 2/3. Note that P(é, 0) = 1. 


It can easily be verified by direct differentiation that if P is a solution of Eq. (16), a 
more general solution can be obtained by setting 

6, £) = | K(é — 7)P(r, £) dz, 
where K is an arbitrary function. In the example considered here we require that 


6,(&, 0) = 4n(2 — &) + 0, ; 


thus we set 
K(é — r) = —4, 
a constant, and @, is given by 
0(€, t) = (0. + 8»)PE, t) — 4v | P(r, $) dr. 


7™R. V. Churchill, Modern operational mathematics in engineering, McGraw-Hill Co., New York, 


\. 


1944, p. 299. 
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In terms of 7 then the first approximation to the temperature distribution in this ex- 


ample has the form 


O(E, n) = (0,, + 8v){PlE, 2v(1 + m)] + Plé, 2x(1 — »)]} 


E 
_ [ {P[r, 201 + »)] + Plr, 2x1 — 2)]} dr 


2 


v-2 
2 





— Sve -2v(1+m) + ” Salta + 4én? . ile 
= 


‘ ~1/3 
where 2p = ¢ 


THE THICKNESS OF A SHOCK WAVE IN AIR* 
By A. E. PUCKETT anp H. J. STEWART (California Institute of Technology) 


1. If the viscosity and heat conductivity of a gas are neglected, it is easy to show 
that a one-dimensional supersonic flow of this gas may be interrupted by a discontinuity, 
or shock wave, across which the velocity jumps to a subsonic value (see [1] and [2]f). 
In a real gas, viscosity and heat conduction may be negligible through a large part 
of a supersonic flow field but clearly must become important in the neighborhood of 
the large velocity gradients associated with a shock wave. If these effects are considered 
in the equations of motion, it is seen that the deceleration corresponding to the velocity 
jump through a shock wave must actually occur in a finite, although short, distance. 

The first theory of the thickness of a shock wave was given by Rayleigh [3] who 
assumed that the fluid was a thermally and calorically perfect gas, that the viscosity 
coefficient « was constant and that the heat transfer coefficient \ was zero. Taylor [1] 
considered the effect of a constant, non-zero value of \ by assuming that the velocity 
jump across the shock wave was small compared with the local speed of sound. Becker 
[4] noticed that a solution could be obtained without linearization in the special case that 
the Prandtl number ¢ = Cpyu/d = 0.75, where C> is the specific heat at constant pressure. 
In all of these investigations the physical constants, Cp , » and A, were considered to 
be constant, and the linear (Navier-Stokes) theories of the viscous stress tensor and the 
heat flux vector were used. All of these investigations showed the width of a shock wave 
in air at normal conditions to be extremely small; for strong shock waves the width 
was computed to be less than the molecular mean free path. In view of these results, 
several writers (see [5] and [6]) have discussed the influence of deviations from thermo- 
dynamic equilibrium and of the higher order, non-linear (Burnett) terms in the viscous 
stress tensor and heat flux vector. 

In the present note the problem of the width of a shock wave in air is re-examined 
using a thermally perfect gas and the linear (Navier-Stokes) theories of the viscous 
stress tensor and the heat flux vector. It is found that the solution can be carried 


*Received Dec. 13, 1948. 
+Numbers in the brackets refer to the bibliography at the end of the paper. 
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through with the physical parameters, Cp , u and X, given as functions of the absolute 
temperature 7', provided o = 0.75. Since the experimental determination of the value 
of \ for air at high temperatures is very difficult, accurate data on the variation of ¢ 
do not exist; however the available data show that o is nearly constant and varies from 
about 0.77 at 0° C to 0.72 at 300° C. The special solution considered here is thus a reason- 
able approximation. If the variation of the physical parameters is taken into account, 
the shock wave width, even for the limiting case of an infinitely strong shock wave, is 
found to be several molecular mean free paths; consequently the importance of the 
Burnett stresses and of the deviations from thermodynamic equilibrium is greatly 
reduced, although not eliminated. 

2. It is assumed that the gas obeys the perfect gas law, P = pRT, where P and p 
are the pressure and density, respectively, and F# is the gas constant. This implies that 
the enthalpy h is a function of the temperature only and that Cp = dh/dT. The con- 


tinuity equation is 


Op 0 
- L , = ( ) I l 
at * ax, (1) 


where u, is the velocity vector. The equations of motion are 


Ou; Ou; 1 oP 1 oO 
= + U; > SS ee i F (7:5); (2) 
ot Ox; p Ox; p Ox; 
where the stress tensor 7,; is given by the linear theory as 
Ou Ou 2 OU, 
ra = wa + 3) — 5 obs 2. (3 
; »( 2 Ox, _— OX, ) 


The corresponding equation for the conservation of energy is most easily derived (see 
[7] or [8]) in the form 
) <i 0 ; D ( P UU fe) , 
d + —(7;,u;) = h—-—-+—] + — (Pu). (4) 
F ( ( Ox; (rast p Dt p 2 Ox; 


Ox Ox 


A simpler form may be obtained by using the continuity equation to transform the 


last term of Eq. (4). It is thus seen that 


re) ar) 0 —_— D ( ut) oP si 
ax, (n az.) * ax, Tit) = 9 Dy ht" at’ 9) 


If these general equations are applied to a plane shock wave, the motion is reduced 
to a steady state motion in a coordinate system stationary with respect to the wave. 
The only velocity component, u, is normal to the wave front; so the flow is one-dimen- 
sional in x, the distance normal to the shock wave front. Under these conditions the 
continuity equation can be integrated to 

pu = m, (6) 


where m is the mass flow constant. The only non-vanishing stress component is 7,, 
(4u/3) (du/dx). Consequently, the equation of motion is 


d »_ 4 au) = 7 
dx (mu 8 37° dx] ~ 8. (7) 
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This can be integrated to 


4 du 


g ha = mu + P — mB. (8) 


Since the velocity gradient must vanish asymptotically far ahead of the shock wave 


(section 1) and far behind the wave (section 2), 


mB = p\ui + Pi = poz + Pe. (9) 
The energy equation (5) becomes 
d (4242) = 2 (yl 4 4 yu tt 
m = (i, + 9 ) a d de + 3 Hua). (10) 
Since h h(T) and Cp = dh/dT, 
dT dh p dh 
“dx  Cpdx odx (11) 


If the Prandtl number is constant, Eq. (10) can thus be written 
| 


d a pee} xd | at ‘)] 
” dx (i 7 9" ~ dx E dx h 4 3 ou : (12) 


This can be integrated at once to 


9) 
h+ = w= At z (i += ou) (13) 


om dx 3 


| 
| 


Since the gradients must vanish at sections 1 and 2, 
A=h+1/2uW =h,+1/2u3. (14) 
Furthermore, if ¢ = 0.75, Eq. (14) shows that the appropriate integral of Eq. (13) is 
h+1/2w=A (15) 
The simple shock theory is obtained by solving Eqs. (9) and (14) together with 
the equation of state in the form pu = mRT and the enthalpy relation h = h(T). 
The theory of the shock wave thickness is obtained by solving Eqs. (8) and (15), 
together with pu = mRT and h = h(T). An interesting feature of this solution is ob- 
tained by considering the point at which the viscous stress [Eq. (8)] is a maximum. At 


this point 


0 = m £ (u + RT). (16) 
ax u 
By Eq. (15) the condition of Eq. (16) becomes (since Cp = Cy + R) 
u? = CERT (17) 
Ce ’ 


i.e., the maximum viscous stress occurs at the point where the speed is equal to the 


local isentropic speed of sound. The fluid properties at this condition will be denoted 


by a re. , ete. 
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The solution of the velocity distribution through the shock wave can be written 

formally, after substituting pu = mRT in Eq. (8), as 

3 o / udu _ 

qr * i + RT — Bw’ (18) 
where the enthalpy, and thus 7 and uy, is known as a function of u by Eq. (15). In 
general, numerical procedures must be used to integrate Eq. (18). It is of interest to 
note that a solution could also be obtained for any value of a; the coordinate x could 
easily be eliminated between Eqs. (8) and (10), leaving a first order, first degree differ- 
ential equation for 7’(u) which could then be integrated numerically as the replacement 
for Eq. (15). 

3. As an example of this theory, consider the case in which the gas is calorically 
perfect, so that Cp is constant, and in which uz is also constant (Becker’s problem). 
For this case Eq. (15) is 

CpT + 1/2u’ = A, (19) 
and Eq. (8) can be written 
4u du i: + ‘) (2 = t) 
b. y.._ , _ | a py aoe 
i lll Dy Bu + ; ; (20) 
where y = Cp/Cy . By the relations of the corresponding simple shock wave theory, 
UU, = uy, B= (uy, + u)(y + 1)/2 and A = wu(y + 1)/2(7 — 1). Thus 


— —. oe _ 
is + De Soe (u u,)(u U2). (21) 


The integral of Eq. (21) is 
om “8 E log ——"" — log ~—* He |, (22) 


Bly + 1)m (uy, — ue) Lue Ux — Uy Ux — Ur 


where z, = 0. Let u = u/u, 80 Uju, = 1; then 
\ 
R, = pu? = —— a — (i log (= — ") — log (= — ')\ (23) 
Ke 3(ty + YD — 1) u, — 1 u, — 1 


where R, is a sort of a shock wave Reynolds number. Since the dimensionless approach 
velocity U, is a function of the approach Mach number, this expression defines the curves 
of R, versus u for various values of /, . This relation is shown in Fig. 1. It is of interest 
to note that this solution for the velocity distribution through a shock wave considering 
both heat transfer and viscosity for ¢ = 0.75 is almost identical with the corresponding 
solution (see [9], p. 651) with no heat transfer (c = ©). For ¢ = the constant factor 
y in the numerator of Eq. (23) is deleted; thus, the effect of heat transfer increases the 
shock wave thickness by the ratio y : 1. 

The shock wave thickness resulting from this computation is, of course, indefinite 
in the sense that % approaches %, and % asymptotically on the two sides of the shock 
wave. However an arbitrary measure of the thickness of the region within which the 
viscosity effects are important may be defined as 


U; pinata U2 _ 
a (—du/dx) (24) 
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For the case that » is constant, (du/dx), is the greatest negative velocity gradient; if 
uw is not constant, it is the velocity gradient at the point of maximum axial viscous 





stress. Even in this latter case Eq. (24) is still a reasonable index of shock wave thickness. 
By Eq. (21) the above equation becomes 
plz d - 8y (“ + us) (25) 
bx 3(y + 1) \u, — Us/" 


As an example, consider a weak shock wave in air at normal atmospheric conditions 
with uy, = 1.05 u, , u,/p, = 1/7 em*/sec., y = 1.400 and u, = 30,000 cm/sec. Then 


2.45) 


| 





6 = 3.0 X 10°* cm. This value is very close to Taylor’s estimate [1]. It is also known 
from the simple shock wave theory that @, can never exceed the value for M, =o, 
Urmax = [Cy + 1)(y — 1)]'” = 2.45, approximately. In this case 
Pld on 
———— = 3.70. (26) 
Mx 


This result can easily be interpreted in terms of the kinetic theory of gases. In the 
simple kinetic theory, 1 = (1/2)pcl, where ¢ is the average molecular velocity and 1 is 
the mean free path. For a Maxwellian velocity distribution, ¢ = 0.921 (83RT)'”; so 
for y = 1.40, 

pees ae. (27) 

Mx ls 
The shock wave thickness for %, = 1.05 (4, = 1.05) is thus 43 mean free paths, for 
Ui, = 1.96 (M, = 3) it is 3.2 mean free paths, and even for u, = 2.45 (M, = ~) the thick- 
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ness is 2.5 mean free paths. A more satisfactory measure of the usefulness of this con- 
tinuum theory is the number of molecular collisions during the time of transit through 
the shock wave. Since the mean molecular speed is approximately equal to the speed 
of sound, the number of collisions is approximately equal to the width, in mean free 
paths, divided by the Mach number. Since M = 1 at u = u,, the shock wave thickness 
as computed above may also be interpreted, roughly, as the number of collisions. The 
interpretation given here differs from those given previously in that conditions within 
the shock wave are used to estimate the mean free path. 

The problem of rate of approach to thermodynamic equilibrium in a shock wave 
has been discussed by Bethe and Teller [6]. They conclude that equilibrium for the 
translational and rotational degrees of freedom is attained very rapidly; probably in 
one or a very few molecular collisions. On the other hand, equilibrium with respect to 
vibrational modes and dissociation may require from 20 to many thousand collisions. 
For air, with low stagnation temperatures, such as would exist in a supersonic wind 
tunnel, the energy content of the vibrational degrees of freedom is almost negligible. 
Consequently it appears that the continuum theory for the thermodynamic changes 
in a shock wave should be a good approximation, since it may be in error by at most a 
few mean free paths, or at most a factor of 2 in the estimated thickness of a shock wave 
at high Mach numbers. If the stagnation temperature is high, much larger errors may 
be expected, since the energy associated with molecular vibration and dissociation may 
be important. 

If the effect of variations in the physical parameters is considered, there are two 
sases which are important. The first is the case in which the variation of Cp is neglected 
but the variation of » is considered. This case applies in supersonic wind tunnels since 
the stagnation temperatures are relatively low, being only a small amount above normal 
atmospheric temperatures in most cases. For this case Eqs. (19) and (21) still apply. 
The integrated form, Eq. (23), is no longer correct since u was considered constant in 
the integration; Eq. (25) is still correct, however, provided the value of » at u = u, is 
used. This result was implied in the kinetic theory interpretation. It may be noted that, 
according to the simple kinetic theory used above, the viscosity coefficient varies as 
T’”’. The actual variation is somewhat more. 

If both » and C> are variable with temperature then Eqs. (8) and (24) give 


Pxilxd _ 4 Ui (uy = Ua) ln (28) 


Ux 3 ur - ly bees U(Yx a LRT 
In general, the shock wave equations (9) and (14) must be treated numerically in order 
o determine uw, , u, and 7, for given values of u, and T, . 
to determir , and T’, for giver lu f nd 7 
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REMARKS ON THE MOTION OF ANCHOR CHAINS* 


By R. GRAN OLSSON (Institute of Technology, Trondheim, Norway) 


1. Introduction. The problem of the motion along a smooth, inclined piane of an 
anchor to which the end of a chain is attached (Fig. 1) has been treated in the book of 
S. Timoshenko and D. H. Young [1].t The equation of motion is, however, given in 








Fia. 1. 


such a form that only a relation between the velocity and the displacement is obtained. 
In this paper it will be shown that very simple relations between the time and the 
kinematic quantities (displacement, velocity and acceleration) can be stated by intro- 
ducing the elliptic functions of K. Weierstrass. Very simple expressions for the kinetic 
quantities (momentum, kinetic and potential energy) can also be established. The 
motion is assumed to be frictionless, but it is not difficult to take into account a dissipa- 
tive force, which is either constant or proportional to the moving mass. 

2. Nomenclature. In this paper the same notations as those of Timoshenko and Young 
will be used: a = acceleration of the weight W, [cm sec~’]; B = rectilinear momentum 
[g sec]; C = constant of integration; g = acceleration of gravity [em sec*];m = mass 
of the system in motion [g em™ sec’]; P = force [g]; ¢ = weight per unit length of the 
chain [g em ']; Q = loss of energy by percussion [g cm]; t = time [sec]; 7’ = kinetic 
energy [g cm]; wu = parameter of the elliptic functions; v = velocity [em sec™']; V = 


*Received May 2, 1949. 
t+Numbers in square brackets refer to the bibliography at the end of this note. 
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potential energy [g cm]; W, = weight of the anchor [g]; W = (W. + qx) = weight of 
the anchor and a length z of the chain [g]; x = displacement [em]; x) = W./q = length 
of the chain, corresponding to the weight W, [cm]; a = angle of the inclined plane. 

Further the following abbreviations are used: s = (Wy + qz)/W, 4°79; = 8 4'% = 
(Wo + qx)/Wo ; u = (g/g) dx = mass of the chain element dz [g cm™ sec’]; p(u) = the 
elliptic p-function of Weierstrass; ¢(u) = the elliptic ¢-function of Weierstrass. 

3. Equation of motion and its integration. By using the second law of Newton in 
its original form, 

a (md) = P, 


and introducing the quantities 


l pw ’ 
m = ra o + qx), (1) 


P = Wsina = (Wy + qz) sina, 
we obtain 


ld ‘ 
<< oneP y = ys y 
oa W sin a. (2) 


After multiplying by g dt, we find that 
d(Wv) = Wg sin a dt, (2a) 


corresponding to Eq. (82c) of Timoshenko and Young [1, p. 114]. 
To integrate Eq. (2a), we multiply both sides by Wv and substitute on the right- 
hand side, W, + gx for W and dzx/dt for v. In this way, we obtain: 


Wv d(Wv) = (Wo + qz)’g sin a dz. (2b) 
Integrating Eq. (2b), we find that 
1 ae g , 3s Y ’ 
= (Wv)” = = (Wo + gz)’ sina+ C. (3) 
2 3g 
Let us assume now that when ¢ = 0, the body is at rest at the top of the incline. 
For this initial condition, the constant of integration C in Eq. (3) becomes 


( a 
C=- £. W; sin a, 
3q 
and we obtain 
2 29(W,.+q2)* —We. 
y= asin aa. 4 
i: +2 ”* (4) 
From this equation we can find the velocity v for any position, provided W, and q are 
known. Taking gz as small compared with W, , we find that the expression (4) reduces to 
v = 29x sin a, 


as for a body of constant mass [1, p. 115]. 
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4. Introduction of the elliptic functions of Weierstrass. Substituting 





— = (W. + qx)/Wo, (5) 
de _ Wed 
dt q dt’ 
we obtain from Eq. (4) that 
ay = Pee £8 
c " oF (4a) 


Integration of this equation gives 


ry \1/2 & 
__ 3W) E dé (4b) 





~ (2gqsin a)” J, @ — 1)” 





where & = 1, corresponding to x = 0. 
A further substitution, 


E = 4'%s, (6) 
transforms Eq. (4b) into the following integral: 


_ WH" fede 
~ (gqsin a)” J,, (48° — 1)'” (4¢) 





Following Weierstrass, we introduce the elliptic integral [2] 


‘i ds - 
U ion! / (48° oe 982 ms er aa ? (7) 





where g, and g; are called the invariants of the doubly periodic function (the so-called 
p-function) 

s = p(u), (7a) 
»(u) being the inverse function of the elliptic integral (7). By means of (7a) we can 
express the time ¢ in the following way: 


6W, vale! ‘ad o eo 
desis = 2. a)? (] p(u) du — / p(u) an}; (4d) 








further, by introducing the ¢-function of Weierstrass defined by [2, p. 183} 
tu) = — | p(u) du, (7b) 


we obtain (7%) = W,./q) 
3.086 (2x\'” 
"Gane (*) [s(@u) — F(%)]. (4e) 
The connection between the parameter u and the displacement z is given by Eqs. (5), 
(6) and (7a): 
— Wo 


. [4'"p(u) — 1] = 2, [4° p(u) — 1). (5a) 


x 
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Thus we have a description of the time ¢ and the displacement x given by the same 


parameter w. 
The square of the velocity can now be calculated from Eq. (4): 


2 : ’ ; ( wo 1 Atl VY, : , 2 
y = - Wy, sina E(u) == gh —sin a p’(u) - (4f) 
3q 4°°°y(u) 6q p(u) 
the velocity itself is 
4/6 - a2 YP’ - raya p’( 
v= B72 (gx sin a)'”” re = 0.5143 (gay sin a)'”” a ; (4g) 


The acceleration can be calculated from Eq. (2), which after differentiating with respect 


to t, yields 


s v+W “ = gW sina. (8) 
Recognizing that 
dw dx? iy 
we obtain the acceleration from Eq. (8): 
3 , 
a= “ = gsina — _ v = gsin of —- wet a ‘ sin a vid. (8a) 


We observe that [see Eqs. (5), (6) and (7a)] 


Wo + az = Wot = W.4'”s = Wot’ p(u); 


x 
NN % 
JO | 
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hence the final expression for the acceleration is 


ee i 1 p(w)’ = -3 
a = gsin of 6 sa = g 5mm a{2 + pu). (8b) 


In the last expression we have used the relation between the p-function and its de- 
rivative p’(u)[2, p. 182], viz., 


p’(u)? = 4p(u)*® — 1. (9) 
Comparing Eqs. (4b) and (7) we see that in the present problem the invariant g2 
is equal to zero; for this case the p- and ¢-functions are tabulated by A. G. Greenhill [3]. 


Thus we can calculate the diagrams for displacement, velocity and acceleration depending 
on time ¢; they are presented in Figs. 2, 3 and 4, respectively. It can easily be proved 


VY 
negara 








AT 
% 


7 2 3 
Fia. 3. 


that the curve in Fig. 4 has a horizontal tangent at ¢ = 0 and a point of inflection for 
p(u) = (5/14)'” = 0.710, corresponding to u = 1.275. 

There is an asymptotic value of acceleration equal to g/3, and this is obtained prac- 
tically with u approximately equal to 0.5. The constant value of acceleration corre- 
sponds to a straight line for velocity and a parabolic curve for displacement. 

5. Presentation of the kinetic quantities. Besides the kinematic quantities, it is 
interesting to present the kinetic quantities—momentum and energy. 
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For the momentum we have: 


_. © head Wo + qa dx Wo jis oi cS; «tai 
B= : y= : ale 4 gi72 (go sin a) ’~p’(u), 


Z 1/2 
B = 0.816617(2sin a) p’(u); 
( 


a 


| ba S/AIA 


Lo | 








4 
Se - 
] 2_ 
|S: S77 

t , > 

¥ 2 3 - — 

Fia. 4 
and for the kinetic energy: 
= W . Wo qx aye Cu)” We p’(u)” 
T=—y =- ot 92 { W, sina 2, = —— sing — . 
29 12q p(w) 6 q p(w) 


Introducing p’(w) expressed by p(w) (Eq. (9)), we obtain 


Wo - a 
f = O:230 sin af4p(u)* — p(u)~"], 
q 


and for large values of the parameter uw, 


T = 0.840W 2 sin ap(u)’. 


Finally, the potential energy 


Vy = — (Wer a qx z) sin a= ~2 sin a W, + 42), 


Inserting the expression for the displacement x (Eq. (5a)), we find that 


V=- = sin al2 Q' pu — De 


(10) 


(11) 


(lla) 
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6. The transformation of energy of the system. It can be shown that the theorem 
of conservation of energy in the usual form does not hold in the present problem. To 
the energy equation we must add a term, involving the loss of energy by percussion of 
the links of the chain. The calculation has the following form. 

The potential energy is obtained by integration, 


V = Wosina(z, — 2) + 3 sin a(xi — 2°), (13) 


and, as above, the kinetic energy is 


rp me _ ley 2 
T= 7 = 2g (Wo + qu)’. (14) 
Differentiating with respect to time yields 
dV ’ _ a 
=a —(W o + ga) sin a 7 ‘ (13a) 
dT v " dv . 
as | 2c o + qa) di + g@ | (14a) 
and after addition, 
dV , aT : ; v , dv : 
wipe = _ ys ) ay 0 “ae ) ° 5 
aT as (Wo + gx) sin a + 3, | 2c + qx) ut @ | (15) 
Equation (2) may be written as 
l ow yw. 92 _ mp li 
; (Wo + qx) rT + “. (W, + qx) sina = 0. (2c) 
If we add the term 
aq _ 43 . 
dt ~ 2g‘ (16) 


to Eq. (15), we obtain 
dV , aT , dQQ_ F . i. , | ‘ 
qi + a +  * - (Wo + qu) 7 + “iy (Wo + qx) sina |, (15a) 
where the terms in the brackets correspond to Eq. (2c). Thus Eq. (15a) may be written 
as 


dV , aT , dQ _ . 
dt dt dt (15b) 


or, after integration, 


V+T7+ Q = const. (17) 
The quantity dQ may be written as 


dQ = 41? dy (16a) 
29 
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and interpreted as the loss of energy by non-elastic percussion of the links of the chain. 


If we write 
p= 4 dx, 
g 
where yu indicates the mass of the chain element dz and 


’ 


dQ = 50" (16b) 


u is also called the reduced mass by percussion of the chain elements (Carnot). 


REFERENCES 
(1) S. Timoshenko and D. H. Young, Advanced dynamics, McGraw-Hill Book Co., Inc., New 


York, 1948, pp. 114-116. 
(2) P. Frank and R. v. Mises, Die Differential- und Integralgleichungen der Mechanik und Physik, 


2. ed., M. S. Rosenberg, New York, 1943, p. 180. 

(3) A. G. Greenhill, Proc. Roy. Artillery Inst. 17, 181-226 (1889); or E. Jahnke and F. Emde, 
Tables of functions, Dover Publications, New York, 1945, pp. 102-104. (The value of ¢(u) = 3.7174 
for r = 31° is incorrect. A better value is ¢ = 3.7952.) 


A NOTE ON SUPERSONIC FLOW IN THE TREFFTZ PLANE* 
By JOHN W. MILES (University of California, Los Angeles) 


Introduction—The general problem of linearized, supersonic flow has been treated 
from the standpoint of Fourier transforms in an earlier paper.’ The methods set forth 
therein will be ysed to calculate the flow conditions far downstream from a lifting sur- 
face. If x is the streamwise coordinate, cf. Fig. 1, the downstream plane x = = is desig- 
nated as the Trefftz plane.” 

This same problem has been treated by Lagerstrom and Graham,’ who took ad- 
vantage of the fact that the flow in the Trefftz plane satisfies Laplace’s equation. These 
same authors also made calculations in regions not infinitely removed from the wing 
by making use of conical flow methods. 

Statement of problem—The vector perturbation velocity due to the presence of the 


lifting surface is given by 


q(z, y,z) = UV o(z, y, 2), (1) 
where ¢ satisfies the linearized equation 
Pzz = Puy + Piz) (2) 


*Received Feb. 11, 1949. 

1J. W. Miles, On linearized supersonic airfoil theory, No. Amer. Avia. Rep. AL-801, Nov., 1948. 

2E. Trefitz, Zur Prandtlschen Tragflaechentheorie, Math. Ann. 82, 306-319 (1921). 

3P, A. Lagerstrom and M. E. Graham, Downwash and sidewash induced by three dimensional lifting 
wings in supersonic flow, Douglas Aircraft Rep. SM-13007, April, 1947. 
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the Mach number being taken as 2'””. The pressure coefficient is given by 


esi, y) = PAS th OM) — BE OF) 8 acme — 1)(e, wd, 3) 
v(x, y) = $2(2, y, 0+). (4) 


Only points in the upper half-space (2 > 0) are considered, by virtue of the (assumed) 
asymmetrical nature of the flow. The problem to be solved is the explicit specification 


Z 








Fia. 1. Coordinate system and lifting surface in the (x,y)-plane 


of the velocity potential ¢(, y, z) from a knowledge of y(z, y) over a lifting surface 
o lying in the (z,y)-plane* (ef. Fig. 1). 
Solution—A solution to Eq. (2), which reduces to y(z, y) at z = 0+, by virtue of 
the Fourier transform relations, and vanishes identically upstream, is given by 
l . e ; Z — , 
$(z, y, 2) = 73 | du / dp [ dé [ dn{(iu)~* + 25(u)]y(E, n) exp {i[u(x — &) 
7 Ja J 2 J—@ J —a (5) 


+ oy — n) — (u’ — v’)'e)}, 


where the path of integration is indented over and under the branch points + | | 
and — |u|, respectively, and under the pole » = 0. The presence of the term 7é(u) 
(where 6(u) is the Dirac delta function) in the kernel is required in order that ¢ will 
vanish at x*= —o. The terms (iu)~* and 76(u), on the other hand, give equal contri- 
butions to the potential at z = +, the asymptotic behavior of the solution at x = © 
being determined by the behavior of its Fourier transform at « = 0 (since the integral 
in the complex y-plane can be closed in» = + i). Thus, it is found that 


(o,y2=5-[ & | a] dmG)expliy— lela. © 


‘y(x, y) of course vanishes for (zx, y) outside of o, since, by hypothesis, only o can support a pressure 


discontinuity. 
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While, in many cases, it may be expedient to carry out the integrations with respect 
to £, » before integrating with respect to v, carrying out the latter integration yields” 

2 @ a y(é, n) 
¢(o, y,z) =- [ dé | dy =. (7 
( ‘ tlc ts FFE 28 ) 


This result may also be obtained (or interpreted) through the method of source distri- 


butions. In the plane z = 0, Eq. (7) reduces to 


a @ 


| 


v—@ 


Ry 


lim ¢(, y,z) = 5s [ dny(é, n)d(y — ») = [ déy(é, y). (8) 

It is of interest to remark that the results, Eqs. (6)-(8), are independent of Mach 
number and are, therefore, valid for all Mach numbers for which linearized theory is 
admissible. In view of Jones’ results for low aspect ratio airfoils,” this state of affairs is 
not unexpected. 

5B. O. Peirce, A short table of integrals, Ginn and Co., New York, 1929, p. 64. 

6R. T. Jones, Properties of low-aspect-ratio pointed wings at speeds below and above the speed of sound, 


NACA TN 1032, March, 1946. 


A THEOREM CONCERNING THE POSITIONS OF EQUILIBRIUM 
OF A RIGID BODY* 


By L. A. MacCOLL (Bell Telephone Laboratories) 


The following interesting theorem is an immediate consequence of three theorems 
which are well known in certain mathematical circles. However, since, as far as the 
writer is aware, the theorem has not previously been stated explicitly, it seems worth 
while to present it in this note. 

Theorem. Let B be a rigid body, of any form, which is free to rotate about a fixed point P. 
Let B be subjected to forces which are derived from a potential energy function V, which is 
an arbitrary single-valued function, of class C*, of the coordinates of B. Then B has at least 
four distinct positions of equilibrium. 

We first note that the positions (or configurations) of equilibrium of B are the points 
in the configuration space at which the function V, which is defined over the configura- 
tion space, is stationary. Then the theorem is an immediate consequence of the following 
three theorems. 

1. The configuration space of a rigid body which ‘s free to rotate about a fixed point 
is homeomorphic with real three-dimensional projective space. 

2. The number of distinct stationary points of a real single-valued function, of class 
C’, defined over a closed manifold M is at least equal to a certain topological invariant of 
M, called the category of M. 

3. The category of real n-dimensional projective space is n + 1. 

The first of these theorems, in one or another of various equivalent forms, is widely 
known. In effect, it is given as an exercise on page 56 ef Seifert and Threlfall’s Lehrbuch 
der Topologie. 


*Received Aug. 8, 1949. 
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The second theorem is stated and proved on page 91 of Seifert and Threlfall’s Varia- 
tionsrechnung 1m Grossen. 

The third theorem is stated and proved on pages 42 and 43 of Lusternik and Schnirel- 
mann’s monograph Méthodes Topologiques dans les problémes variationnels, Premitre 
Partie, Paris, Hermann et C**, 1934. (Actualités Scientifiques et Industrielles, No. 188.) 


NOTE ON THE KERNEL exp(— | x — y| )* 
By HARRY POLLARD (Cornell University) 
In a recent issue of this Quarterly a theorem was stated which, after some trivial 
changes and correction of a typographical error, reads as follows." In order that the func- 


tion f(x) have the form 


, 1 ee —Iz ! ‘ 
(1) fas e *"g(y) dy, z> 0, 


“0 
it is necessary and sufficient that 


(a) f(0) = f'(0); 
(b) g(x) = f(x) — f(x); 


g(x) = O(e), x >, for some d < 1. 


Now the example f(x) = e* shows that the conditions are certainly not sufficient even 
if (c) is replaced by the condition 
(c’) the integral in (1) exists, 
which the author regards as equivalent to (c). 

On the other hand the example g(x) = e’(x” + 1)~' shows that the condition (c) is not 
necessary either. It also shows that (c) and (e’) are far from equivalent. 

The following is a correct version of the theorem. The only difficulty is the discovery 
of a condition to replace (c). Once this is done the proof is straightforward, and is there- 
fore omitted. 

THEOREM. Jn order that f(x) have the form (1) where g(y) is a prescribed function 
integrable on each finite interval it is necessary and sufficient that f’(x) exists and is absolutely 


continuous on finite intervals, and moreover that 


(a) f(0) = f'(0); 


(b) g(x) = f(x) — f’’(x) for almost all x > 0; 
(c) f'(x) = o(e*), r > &. 


z 


temarks. (i) The condition (¢) cannot be replaced by O(e*) as the example f(x) = e 
shows. (ii) If the word “prescribed”’ and condition (b) are omitted simultaneously the 
theorem remains true. 

*Received August 12, 1949. 

1H, P. Thielman, Q. Appl. Math. 6, 443-448 (1949), Theorem I. 








174 BOOK REVIEWS [Vol. VII, No. 4 
BOOK REVIEWS 


Integraltafel. Erster Teil: Unbestimmte integrale. By Wolfgang Grébner and Nikolaus 
Hofreiter. Springer-Verlag, Vienna, 1949. viii + 166 pp. $5.40. 


This table of indefinite integrals appears to be more complete than others generally available. There 
are three sections of the book entitled Rational Integrands (21 pages), Algebraic Irrational Integrands 
(85 pages), and Transcendental Integrands (60 pages). 

Algebraic Irrational Integrands contains integrals having integrands as follows: 
Rational functions of x and (ax + b)'/"; 2 and (ax + b)!/*; 2 and [(axr + b)/(ex + d)}"/"; 2, (ax + b)'/2, 
rand (ax? + 2br + c)'/? with a number of special cases of this last type. Irrational inte- 


The section on 


and (er + d)! 
grands which can be transformed into rational integrands, Elliptic integrals in Legendre canonical form 
and in Weierstrass canonical form. Rational integral functions of z and y = (aor? + 3a,x2? + 3agr + a;)"/2, 


Reduction to the Legendre canonical form. Rational integral functions of z and y = (aer* + 4ayr? + 


+ a,)'/?, Reduction to the Legendre canonical form. Rational integral functions of x and 


6a | 1a r 
a;)]'/3, Reduction to the Weierstrass and 


y = (aor? + 3a,r? + Baer + a;z)'/? = [ag(x ai)(x — ae)(x 
Legendre canonical form. Rational integral functions of x and y = (2? + 1)'/3, Reduction to the Legendre 
canonical form. Hyperelliptie integrals. 

The section of the book on transcendental integrands contains integrals having integrands of the form 
R(exp Ar), f(x) exp Az, f(x) exp (ax? + 2br +c), S(log x), R(x) log . f(x) log "Zr g(x), R (sin x, cos x P 


R(sin (ax + b), cos (ex + d), +++), 2? sin ™x cos "2, exp ax sin ™br cos "cr, R(x, exp az, sin br, cos cr), 
! 4 i ‘ sit tar oa “* 
R(S'" (az? + 2ber +c), x), R(x, are V1" x), Ria, arc @Y x), R(Sinh xz, Cosh x), R(Sinh (ar + b), 
iia ‘ +: . . Sink 
d), +++), 2? Sinh ™r Cosh "z, R(Sinh (az + 5b), sin (cr + d), +--+), R(x, are Conk x) and 


Cosh (ex 4 
R(z, are fare r). Finally there are integrals of Weierstrass elliptic functions and Jacobi elliptic functions, 
This table of integrals is labeled Part I of two parts. Part two on definite integrals is said to be in 


preparation. 
Roun TRUELL 


Fluid dynamics. By Victor L. Streeter. McGraw-Hill Book Company, Inc., New York, 
Toronto, London, 1948. xi + 263 pp. $5.00. 


According to the preface “this book is planned to introduce the reader to the general theory of fluid 
flow” and is expected to be useful as a text for a second course for seniors or beginning graduate students. 
The mathematical background assumed extends through elementary calculus, additional mathematics 
being developed as needed. The author states he has made every effort “to clarify the concepts and to 
include those exasperating steps in deviations which are usually omitted.” 

Except for a part on viscous flow, the book is restricted to the study of an “ideal” fluid, by which the 
author means incompressible as well as nonviscous. Although some of the equations in the text hold for a 
compressible fluid, this subject is essentially not touched. No wave motion is discussed and turbulence 
receives but a brief mention. No mention is made of the Froude number or of the effects of surface tension. 
Thin airfoil theory and lifting line theory are not presented. Hence the title and preface are completely 
misleading as to the real scope of the book. Within the field of incompressible nonviscous flow the treat- 
ment is quite complete, particularly in two-dimensional flow where complex variable theory is used 
extensively. Flows described by inverse transformations and jet and wake flows are well presented. 

The conceptual presentation is basically limited by the scope of the book. Since no thermodynamics 
is used the concept of the energy equation as distinguished from Bernoulli’s equation is omitted. Simi- 
larly the basic concepts of barotropy and piezotropy and their relation to circulation are not included. 
The fundamental fluid dynamic concepts required in compressible flow, dynamic meteorology, chemical 
or heat power engineering, and open channel flow are all absent. 

Within the scope of the book the concepts are well treated as a rule. Weaknesses occur, however, 
generally as omissions. Uniqueness theorems are given for cases where the velocity potential is single- 
valued with no mention of the important cases where the potential is multi-valued. Euler’s momentum 
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integral is omitted. In the theory of vortex motion the law for the motion of an individual vortex is not 
given, and hence there is no discussion of the effect of solid boundaries. The author discusses qualitatively 
the motion of one or two circular vortices apparently oblivious of the fact that the self-induced velocity 
of such a vortex is infinite. The logic used in the discussion of percolation is faulty. 

(ll in all, this book is a good text for the potential flow of an incompressible fluid with an adequate 
first treatment of viscous effects. The mathematics is restricted in level to the complex variables developed 
in the text. The short cut notation of vector analysis is not used and the “exasperating steps” are in- 
cluded. As such it should serve as a senior text. The plan of the book to introduce the general theory of 
fluid flow is not fulfilled, however, and the title should have been worded differently. 


WaLtace D. Hayes 


Tables of inverse hyperbolic functions. By the Staff of the Computation Laboratory. 
Harvard University Press, Cambridge (Mass.), 1949. xx + 290 pp. $10.00. 


The volume contains the following four tables: I) tanh! z for z = [0(.001) 0.5(.0005) 0.75(.0002) 
0.90(.0001)0.95(.00005) 0.975 (.00002) 0.99(.00001) 1], II) sinh! z for x = [0(.002) 3(.005) 3.5], ITT) 
cosh~! x for z ([1(.00001) 1.004(.00002) 1.01(.00005) 1.025(.0001) 1.05(.0002)1.15(.0005) 1.4(.001) 
1.8(.002) 3.5], and IV) sinh~! x and cosh! z for z = [3.5(.005) 6.5(.01) 15(.02) 35(.05) 65(.1) 150(.2) 
350(.5) 650(1) 1,500(2) 3,500(5) 6,500(10) 15,000(20) 23,000]. The first and second differences are given 
in all these tables; the functions are tabulated to 9 decimal places. 

W. PRAGER 


Tables of Bessel functions of fractional order. Prepared by the Computation Laboratory 
of the National Applied Mathematics Laboratories, National Bureau of Standards. 
Volume II. Columbia University Press, New York, 1949. xviii + 365 pp. $10.00. 


The first part (pp. 1-275) contains tables for J,(x) for »y = +1/4, +1/3, +2/3, 43/4 and z ranging 
from 0 to 13 for the negative values of vy and from 0 to 25 for the positive values of v. For small values 
of up to 1.00 for » = —3/4, —2/3, to 0.80 for »y = —1/3, —1/4, to 0.60 for »y = 1/4, 1/3 and to 
0.50 for »v = 2/3, 3/4) the interval is 0.001; for larger values of z the interval is 0.01. The second part 
(pp. 277-332) gives the values of e~7/,(r) for the same values of v and for x = [25(.1) 50(1) 500(10) 
5,000(100) 10,000(200) 30,000]. All functions are tabulated to ten decimal places or ten significant figures. 


A number of auxiliary tables facilitate interpolation in the z and » directions. 


W. PRAGER 


Electromechanical transducers and wave filters. By Warren P. Mason. Second Edition. 
D. Van Nostrand Company, Inc., Toronto, New York, London. xii + 419 pp. $6.00. 


The second edition of this book differs from the first mainly in the addition of new topics in which 
interest has grown as a result of many developments during the war. Those who have read the first 
edition will recall that the subject of linear vibratory systems was developed from the point of view of 
electrical network theory. After a review of electrical network theory the governing equations for lumped 
mechanical, acoustic, elastic, piezoelectric, and miscellaneous electromechanical systems are discussed 
and the electrical network equivalences are developed. 

In the present edition the new material is added in the form of a series of appendices, with the 
paragraphs numbered in such a way as to follow in a natural way some section of the first edition. The 
new material covers some one hundred pages and treats in a more general manner than did the first 
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edition, the motion of a vibrating bar, wave propagation in a viscous medium, and piezoelectric phe- 
nomena. In addition, there is discussion of the electromechanical properties of systems incorporating 
gyroscopic elements, a discussion of the behavior of electrical wave guides and acoustic delay lines as 
both lumped circuit elements and radiation devices. 

The text is extremely clear and is quite readable for engineers, mathematicians and physicists. 
It is not only useful as a discussion of this subject, but also equally valuable as a concise reference on 
electromechanical systems and electrical-mechanical analogies. 


J. A. KRUMHANSL 


Advanced dynamics. By 8. Timoshenko and D. H. Young. Volume I. McGraw-Hill Book 
Company, New York, Toronto, London, 1948. xii + 400 pp. $5.50. 


This book deals with methods of solution of dynamical problems in engineering. Chapter I gives 
a detailed presentation (106 pp.) of various methods of graphical and numerical solution of certain 
non-linear problems, which can be treated in terms of a single particle. Chapter II (84 pp.) deals with 
certain problems of systems of particles, and treats particularly the theory of balancing of reciprocating 
engines. In both of these chapters only elementary notions of single particle dynamics, or of energy 
and momentum relations of a system of particles are employed. In Chapter III (54 pp.) generalized 
coordinates are introduced, and Lagrange’s equations are derived. A number of problems of systems 
involving constraints are treated, and an approximate solution of a problem of forced vibration of a 
system with a non-linear spring characteristic is based on Hamilton’s principle. The fourth chapter 
(82 pp.) introduces the general theory of small vibrations of a system of particles, with detailed dis- 
cussion of vibrations of systems with two degrees of freedom. There is a short discussion of stability 
of vibrations about a state of steady motion, and an introduction to the theory of variable-speed vibra- 
tion dampers. Finally, in the fifth chapter (55 pp.) Euler’s equations are developed and a treatment 
of gyroscopes is based on them. Although the content of this book is not the “advanced dynamics” 
that a mathematician or physicist might expect from the title, the book should prove most valuable 
as a basis for courses in mechanical vibrations and dynamics applied to engineering problems on the 
senior or early graduate level in engineering schools. 

P.S. SyMonps 


Numerical methods of analysis in engineering (successive corrections). A publication 
resulting from a Symposium held at Illinois Institute Of Technology, Chicago, 
Illinois. Arranged and edited by L. E. Grinter. The MacMillan Company, New York, 
1949. xvi + 207 pp. $5.80. 


Most of the chapters of this book originated as papers presented at a Symposium on Numerical 
Methods held at Illinois Institute of Technology in 1948. Since Hardy Cross’s now classical paper on 
“Analysis of Continuous Frames by Distributing Fixed-End Moments” provided an important impetus 
to the application by engineers of numerical methods, this paper is reproduced as Chapter 1, and the 
book is dedicated to Professor Cross. The viewpoint throughout most of the book is practical rather 
than mathematical; the authors are concerned with the quickest methods of obtaining “engineering 
solutions’, and the emphasis is on methods which yield successive corrections. 

Chapter 2 by L. E. Grinter presents a method of solution of plane problems of elasticity by re- 
placement of the continuum by a grid of columns and beams which may be analyzed for example by 
distribution methods. Chapters 3 and 4 by F. S. Shaw and R. V. Southwell, respectively, are concerned 
primarily with relaxation methods applied to boundary value problems; the former reviews techniques, 
and the latter discusses the problem of improving accuracy, contrasting the ‘‘“mathematician’s approach”’ 
by means of more elaborate finite difference formulas with the “engineer’s approach” in which the main 
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concern is to reduce the intervals of the net. In Chapter 5 M. M. Frocht reviews the “linear rosette” 
method for obtaining good initial values in numerical solutions of Laplace’s equation. Chapters 6 and 
7 are concerned with solutions of typical problems of heat conduction by numerical methods; these 
are by L. M. K. Boelter and Myron Tribus, and by G. H. Dusinberre, respectively. 

The concluding group of three papers are concerned with general comparisons between numerical 
methods. Chapter 8 by F. Baron discusses the relations between methods which have been applied to 
problems in diverse fields. In the following chapter N. M. Newmark classifies and compares methods 
which have been applied particularly to problems of elasticity. The final chapter 10, by T. J. Higgins, 
consists essentially of an annotated bibliography of papers in which numerical methods are applied 
to a solution of Saint Venant’s torsion problem; the viewpoint here is somewhat more theoretical, and 
such methods as those of Ritz, Trefftz and Galerkin are briefly summarized as well as the finite difference 
methods with which the preceding chapters were mostly concerned. 

The book provides a valuable cross-section of current applications of certain types of numerical 
methods to engineering problems. Many gaps, of course, exist. The beginner should perhaps be warned 
that numerous misprints are present. 

P. S. SymMonps 


Table for use in the addition of complex numbers. By Jérgen Rybner and K. Steenberg 
Sorensen. Jul. Gjellerups Forlag, K@benhavn, 1948. xiv + 95 pp. $5.50 


The object of the table is to facilitate the calculations of addition or subtraction with complex num- 


bers by tabulating the complex relation 
Re ja =x l 4+ re id 
or Ria = 1 + rio 
The table lists the numerical value of R and the phase angle a of the sum as functions of r and ¢ for 
0 < < 1 at intervals of 0.01 and for 0 < ¢ < 180° at intervals of 1°. a is expressed in degrees and 
decimal fractions thereof. Interpolation differences are given. The authors have prepared this numerical 
table because in a previously published work entitled Nomograms of Complex Hyperbolic Functions by 


Jérgen Rybner it was not possible to construct nomograms with an accuracy sufficient for practical 
calculations. (The tables may be purchased from Scandinavian Book Service, P. O. Box 99, Audubon 
Sta., New York 32, N. Y.) 

Roun TRUELL 


A concise history of mathematics. By Dirk J. Struik. Dover Publications, Inc., New York, 
1948. Vol. I, xviii + 123 pp. Vol. II, 175 pp. $1.50 per Vol. $3.00 per Set. 


This new history of mathematics is a volume in the Dover Series in Mathematics and Physies and 


which has thus far largely included translations into English of volumes in the well-known collection 
known as Sa? ng Géschen. The latter contained a history of mathematics of about the same length as 
Struik’s written by the great historian of mathematies, Heinrich Wieleitner, and published in 1922. It 
was a happy thought to have Struik write a wholly new work rather than simply translate Wieleitner’s, 
because the latter’s Geschichte der Mathematik is “out of date’? and, furthermore, is not that author’s 


best work (his Jathematische Quellenlnicher, for example, are much better). Some of the particular merits 
of Struik’s presentation can best be seen by comparing it to Wieleitner’s. 

One of the first features that strikes the reader of the work under review is the considerable amount 
of bibliographical material, at once leading the student to further reading and showing him that the 
study of the history of mathematics is an active field in which new discoveries are being made, just as 
they are being made in the field of mathematics itself. Struik’s interesting and extended discussion of the 
mathematics of primitive peoples is an important novelty. The section on the ancient Orient, especially 
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rich on the mathematics of the Babylonians, is another important feature. The stress placed by Struik 
on Babylonian mathematics and primitive mathematics (the first treated only summarily by Wieleitner 
and the second ignored by him altogether) is an indication of features of the history of mathematics 
which were not present when Wieleitner wrote his history. The new discoveries in the field of Babylonian 
astronomy and mathematics, due chiefly to O. Neugebauer, have been made largely in the last 25 years, 
and the general interest in anthropology as related to the history of science is also of recent date. Whereas 
Wieleitner divided the modern period by topic, so that one has, for the 18th and the first half of the 19th 
centuries, snippets on each possible mathematical subject, Struik has wisely couched his exposition of 
19th-century mathematics in terms of mathematicians and schools. The result is that Struik’s account 
takes on humane perspectives and is eminently readable, whereas Wieleitner’s history tends to be a 
catalogue. Finally, Struik has included a considerable amount of social history, both ordinary social and 
economic history, as well as an interpretation of the social history of mathematics. This general approach 
places the mathematical history in a context of the greatest interest for the reader, even though at times 
it sets forth a point of view or interpretation with which all readers will not agree. 

Following the opening account of the formation of mathematical concepts among non-mathematical 
people, is a splendid discussion of pre-Greek mathematics, setting the stage for the great accomplishments 
of the Greeks. (Here one must give Struik special credit for his attempt to distinguish between estab- 
lished fact, hypothesis, and tradition.) This chapter is followed by a short discussion of the Islamic period 
which leads him naturally into the revival of mathematics in western Europe, and completing the first 
volume. The second volume, containing three major sections—the 17th century, the 18th century, and 
the 19th century—will probably prove of greater interest to active mathematicians. Nevertheless, in 
contrast to the first volume, in which the major ideas are carefully explained, in this one Struik has 
allowed the mathematics to speak for itself. Only a fully rounded mathematician who knows the major 
fields of mathematics will be able to fully appreciate much of this volume. 

Struik’s book thus clearly fills a long-felt need for a short readable history of the main ideas and major 
individuals in the history of mathematics. It is a tribute to his mastery of the craft of exposition that he 
has been able to include so much related material on astronomy and mathematical physics and so much 
up-to-date information concerning the history of mathematics. Struik has not only read the original 
writings of the great mathematicians whose work he recounts, but he has also mastered the vast secondary 
literature of biography, interpretation, and commentary. 

In a work at once as concise and as comprehensive as this one is, certain omissions are bound to occur 
to each reader. Nevertheless, I feel sure that no one will be able to find a major mathematician omitted. 
Nor in general can anyone but admire Struik’s selection of topics. Anyone who reads this concise history 
of mathematics will, I am sure, join with me in expressing the wish that the work under review be con- 
sidered only an earnest of a major work to come, in which Struik will have the opportunity to explore 


further the social history of mathematics. 
I. BERNARD COHEN 


Numerical calculus. (Approximations, interpolation, finite differences, numerical integration, 
and curve fitting.) By William Edmund Milne. Princeton University Press, Princeton, 
N. J., 1949. x + 393 pp. $3.75. 


Numerical analysis has been sadly neglected in our traditional undergraduate and graduate curricula 
in mathematics because it was not considered worthy of the attention of serious mathematicians. Now 
that the availability of electromechanical and electronic calculators has revived interest in this field, the 
scarcity of adequate textbooks, a consequence of previous neglect, confronts mathematics departments 
which are ready to increase the course offerings in this field. The appearance of the present text is there- 
fore particularly timely. 

The scope of the book is best indicated by the following list of chapter headings: simultaneous linear 
equations, solution of equations by successive approximations, interpolation, numerical differentiation 
and integration, numerical solution of differential equations, finite differences, divided differences, 
reciprocal differences, polynomial approximation by least squares, other approximations by least squares, 
simple difference equations. These eleven chapters fill 348 pages; the remaining pages contain a table 
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comparing the various notations used in the calculus of finite differences, a classified guide to the formulas 
and methods presented in the book, a bibliography, and several brief numerical tables. 

The presentation has been kept elementary; the usual three semester sequence of elementary calculus 
and ordinary differential equations should provide sufficient background for the users of this book. The 
remarkably low price of the cloth-bound volume has been made possible by photo-offset reproduction 
from a carefully typed manuscript. 

W. PRAGER 


Quantum theory of fields. By Gregor Wentzel. Interscience Publishers, Inc., New York 
and London, 1949. ix + 224 pp. $6.00. 


The advanced student of modern physics must certainly welcome this English edition of the most 
complete and probably the best introduction to the quantum theory of fields. The translation from 
Einfiihrung in die Quanten theorie der Wellenfelder (Wien 1943) was made by C. Houtermans and 
J. M. Jauch; it includes a preface to the English edition by G. Wentzel. While the most recent theoretical 
developments in quantum electrodynamics have not been incorporated in the references, the basic 
information is provided to enable the reader to follow the current literature. Some changes have been 
made in the section on the applications of field theory to problems of nuclear physics, and in various 
sections of the book references to recent publications have been added. An Appendix has been added deal- 
ing with the definition of a symmetrical energy-momentum tensor for any kind of a field with an invariant 
Lagrange function. 

The contents of the book fall under the chapter headings: General Principles, Scalar Fields, The Vec- 
tor Meson Field, Quantum Electrodynamics, The Quantization of the Electron Wave Field According 
to the Exclusion Principle, Supplementary Remarks, and the Appendix mentioned. 


Roun TRUELL 


Non-linear problems in mechanics of continua. Proceedings of Symposia in Applied 
Mathematics, Volume I. American Mathematical Society, New York, 1949. vii + 219 


pp. $5.25. 


This volume contains the papers presented at the First Symposium on Applied Mathematics of 
the American Mathematical Society, held at Brown University in the summer of 1947. The subject of 
the Symposium was the same as the title of this book. 

The contents are divided into two groups, one under the subtitle of Hydrodynamics and the other 
under that of Elasticity and Plasticity. The majority of the contributions in the first group are in the 
field of compressible flow, with studies of two-dimensional flows, shock wave problems, and the general 
three-dimensional case. Included also are contributions on incompressible flows with free boundaries, 
turbulence, and boundary layers. In the second group investigations are reported in general elasticity, 


large deflection plate theory, and plasticity. 
W. Hayes 


Tables of the Bessel functions of the first kind of orders fifty-two through sixty-three. By The 
Staff of the Computation Laboratory. Harvard University Press, Cambridge, Massa- 
chusetts, 1949. 544 pp. $8.00. 


The Bessel functions J*.(x), J*°s(x), +++ J®s(x) are tabulated to ten decimal places. The argument 


varies in steps of .01 forO0 < x < 100. 
G. F. CARRIER 





480 BOOK REVIEWS [Vol. VII, No. 4 


Tables of supersonic flow around cones of large yaw. By The Staff of the Computing Section 
Center of Analysis. Under the direction of Zdenek Kopal. Massachusetts Institute 
of Technology, Cambridge, Massachusetts, 1949. xviii + 125 pp. $1.50. 


Some results concerning the flow of a non-viscous gas past a cone with yaw are tabulated in this 
volume. The theory is as yet unpublished, but the compilation of formulas given in the introductory 
pages indicates that this theory utilizes the first two terms in a perturbation expansion in powers of 
the yaw angle. Cones of semi-apex angle 5°, 7.5°, 10°, 12.5°, 15°, 20°, 25°, are considered. The pertinent 
quantities (which define the pressure, velocity, . . . ) are plotted versus the polar angle. The increments 
in this independent variable are .1 near the shock and near the body and are larger in the intermediary 
region. No comments concerning applicability and accuracy can be made until the details of the theory 
are available. 

G. F. CARRIER 


The Aeronautical Quarterly. Published by the Royal Aeronautical Society, 4 Hamilton 
Place, London, W.1. Volume I: Part I (May 1949), 122 pp., 10s. Od., Part II (August 
1949), 72 pp., 10s. Od. 


As is stated in the foreword to the first issue of this new Quarterly, aeronautical science has grown 
more rapidly in recent years than the means of disseminating the knowledge acquired during this growth. 
The new Quarterly was founded as a means of remedying this situation; its Editorial Board encourages 
research workers in aeronautics to submit papers describing new and original work, or papers reviewing 
the progress in some specialized field. 

The first two parts of Volume I contain a diversified group of papers treating, among many others, 
such topics as downwash behind a supersonic wing, flutter of systems with many freedoms, displacements 
of elastic systems under transient loads, and laminar layers in compressible flow. 

The addition of this new Quarterly to the family of journals devoted to problems of engineering 


science is welcomed. 
G. F. CARRIER 


Mathematical foundations of statistical mechanics. By A. I. Kinchin. Dover Publications, 
Inc., New York. viii + 179 pp. $2.95. 


This book forms an excellent introduction to the difficult and important discipline of Statistical 
Mechanics. It is clear, concise and rigorous. There is a very good chapter on the ergodic theorem (with a 
complete proof!) and a traditional but highly lucid chapter on statistical foundations of thermodynamics. 
The treatment is somewhat reminiscent of that of P. Hertz (misspelled Kertz in the book) and is based 
on the notion of the microcanonical ensemble. The scope is extremely limited since the author deals only 
with systems composed of many independent parts (ideal polyatomic gas). It is because of this limitation 
that a variant of the central limit theorem can be successfully applied. One gathers, unfortunately, the 
impression that the author tries to replace once and for all the ‘abstruse’ method of steepest descent by 
the use of the central limit theorem. The reviewer fails to see how the central limit theorem can be helpful 
in cases when non-negligible interactions between the components are present. It is a pity that in a book 
entitled ‘Foundations’ no mention is made of the Boltzmann-Zermelo polemic concerning the meanings 
of the H-theorem (in fact, there is no mention of the H-theorem) and there is no discussion of non- 
equilibrium problems. The translation is very good and only in few places did the reviewer find his 
knowledge of Russian helpful. 

The book should prove useful to teachers of Statistical Mechanics and to mathematicians who 


would like to get acquainted with the problems and methods of this fascinating branch of physics. 


M. Wace 























